Using runs of cells to traverse a ray through a volume

ABSTRACT

Line drawing techniques that employ runs or runs of runs of pixels to draw the line compute line structure information that they use to determine the sequence of runs in the line. This line structure information may be used to compute the positions of a plurality of the runs and then draw the runs in parallel. The line drawing techniques may be also be used with rays in three dimensions. Projections of the ray are made on planes that intersect each other on the ray&#39;s major axis. The line drawing techniques are used to determine cells in the planes that are intersected by the projections. The voxels intersected by the ray are then determined using the cells. Runs of voxels in the ray are used in ray traversals. The volume traversed by the ray is subdivided into encoding runs of voxels that may include one or more significant runs containing voxels whose data will affect the ray. Traversal is done by determining for each run of voxels in the ray whether any of the voxels in the ray run are also in a significant run.

CROSS REFERENCES TO RELATED APPLICATIONS

The present patent application claims priority from the followingapplications:

-   -   U.S. provisional patent application 60/346,741, Peter        Stephenson, Method for visualizing volumetric datasets, filed 7        Jan., 2002;    -   PCT/US02/24711, Stephenson, et al., Methods and apparatus for        determining intersections of a particular line with cells in a        lattice filed 2 Aug. 2002; and    -   PCT/US02/39716, Peter Stephenson, et al., Using line structure        information to enhance line drawing in digital systems, filed 12        Dec. 2002        PCT/US02/24711 claims priority from U.S. provisional patent        application 60/309,926, Stephenson, et al., Process and        apparatus for line drawing, filed 3 Aug. 2001. PCT/US02/39716        claims priority from U.S. provisional patent application        60/341,503, Stephenson, et al., Process and apparatus for        anti-aliased line drawing, filed 13 Dec. 2001, and U.S.        provisional patent application 60/341,194, Stephenson, et al.,        Process and apparatus for line drawing in parallel, also filed        13 Dec. 2001.

In the US national stage, this application will be acontinuation-in-part of both PCT/US02/24711 and PCT/US02/39716. Thisapplication contains the entire Detailed Description of PCT/US02/24711and the portion of PCT/US02/39716 titled Parallel techniques fordetermining the cells of a raster that are intercepted by a line. Thematerial in the Detailed Description, which is new to this application,begins with the section Problems of ray traversal in three dimensions.

BACKGROUND OF THE INVENTION

1. Field of the Invention

The invention relates generally to techniques for determining whichcells of a raster a particular line intersects. The cells of the rasterare represented in the memory of a computer system and the computersystem's processor makes the determination. The cells of the raster maygenerally represent a set of locations. The locations may be pixels in adisplay, and when they are, the techniques may be used to determinewhich pixels in the display represent the particular line, and thus togenerate the line in the display.

2. Description of Related Art

FIGS. 9-11

Systems Using Raster Display Devices: FIG. 9

The flat panel or cathode ray tube display devices typically used withcomputer systems are raster devices, that is, the display area is madeup of a large number of picture elements, or pixels, which are arrangedin a grid. The location of any pixel in the display can be specified byits row and column in the grid. Varying the color and intensity of thepixels in the grid makes the image that the display device displays.

FIG. 9 is a high-level overview of a computer system with a rasterdisplay device. The main components of system 901 are a processor 911with memory 903 to which processor 903 has access and a monitor 915 forwhich processor 911 generates displays. Monitor 915 is a raster displaydevice and as such, has a grid of pixels 917. Within memory 903 arestored bitmap 909, bitmap drawing code 985, and bitmap data 907. Bitmap909 is an area of memory that corresponds to grid of pixels 917. Eachitem of data in bitmap 909 corresponds to a pixel in grid of pixels 917.Drawing code 985 is code for drawing graphical entities such as lines orpolygons in bitmap 909. Bitmap data 907, finally, is data, typicallysupplied by a user program, which bitmap drawing code 985 uses to draw agraphical entity. For example, if a user program wishes to specify aline, it will typically indicate the start and end coordinates of theline in grid of pixels 917 and drawing code 985 will use thatinformation to draw a corresponding line in bitmap 909. Processor 911then reads from bitmap 909 to generate an image for display on grid ofpixels 917. In many cases, a display generator component 913 ofprocessor 911 reads bitmap 909 and produces the actual signals formonitor 915 from the data in bitmap 909. It should be noted here thatthe task of producing a display 917 may be distributed in many differentways across hardware and software components, with the amount ofhardware increasing as performance demands increase.

Representing Lines Using Pixels: FIG. 10

Drawing straight lines is a problem with any raster display device. FIG.10 shows why. FIG. 10 shows a representation in pixels 1001 of the line1003 that is described by the equation $y = {\frac{17}{41}{x.}}$The representation includes those pixels in the grid that areintersected by line 1003. These pixels form a pattern 1004 that istermed the intersection pattern for the line. A line's intersectionpattern depends not only on the line's slope, but also on the locationof its endpoints relative to the grid of pixels. As will be explained inmore detail in the following, the intersection pattern for any straightline has regular features that can be used in drawing the straight linein a raster display or analyzing a straight line that is displayed in araster display.

At its lowest level, the intersection pattern for line 1003 is asequence of pixels. For a given next pixel, there are only twopossibilities: if the current pixel has the coordinates (a,b), the nextpixel has either the coordinates (a+1,b) or (a+1,b+1). Which of the twopossibilities the next pixel has depends on where the line intersectsthe current pixel. To draw a line one pixel at a time, one need onlydetermine for each pixel where the line intersects the current pixel anduse that information to determine which of the two possible positions togive the next pixel.

As is apparent from FIG. 10, the intersection pattern includes groups ofadjacent pixels that have the same y coordinate. Such a group of pixelsis called a run. One such run of three pixels is shown at 1005. Anexamination of the runs in FIG. 10 shows that they have only twolengths: a short length 107, which is here two pixels, and a long length1009, which is here three pixels. The general rule is that the runs ofan intersection pattern will have only two lengths, and these lengthswill be consecutive integers. The lengths of the long and short runs forany given line can be computed as follows:

Within the intersection pattern of the line l:${{l\text{:}y} = {\frac{a}{d}x}},$there are a runs that correspond to the Y-axis size of the lattice. Asthere are only two possible run lengths in a given intersection patternand the possible lengths are consecutive integers we will refer to thelengths as short (s) and long (l). To determine what run lengths arepossible within the intersection pattern, consider that there are dpixels to be distributed among a runs, the distribution being as even aspossible. If we divide up the d pixels into a runs of length$r = \left\lfloor \frac{d}{a} \right\rfloor$we have n≡d mod a pixels remaining, 0≦n<a, which have to be distributedalong the intersection pattern. Therefore in the intersection pattern ofl there are n long runs each with r+1 pixels and a−n short runs with rpixels each.

This can be applied to line 1003 as follows: line 1003 contains 41pixels and 17 runs. The possible run lengths are$r = {\left\lfloor \frac{41}{17} \right\rfloor = 2}$and r+1=3. There are 41 mod 17≡7 long runs 1009 of length three shown inlight gray and 17−7=10 short runs 1007 of length two shown in dark gray.Therefore using a run-based algorithm improves upon pixel-basedalgorithms as only a and not d decisions whether to increase theycoordinate by 1 are necessary.

An examination of intersection pattern 1004 of line 1003 shows that thelong and short runs themselves occur in repeating patterns. Thus, inintersection pattern 1004, there is a repeating pattern of a long run(l) followed by two short runs (s) followed by a long run followed byone short run, or lssls, as shown at 1011 and 1013. In general, thereare four possibilities for the patterns of runs:

-   -   ls⁺, a long run followed by one or more short runs;    -   l⁺s, one or more long runs followed by a short run;    -   sl⁺, a short run followed by one or more long runs; and    -   s⁺l, one or more short runs followed by a long run.

These patterns are termed in the following the shapes of runs. Thus,using this notation, the shape of the runs shown at 1011 and 1013 isls⁺. Moreover, it turns out that the first run in the intersectionpattern of the line l: ${l\text{:}y} = {\frac{a}{d}x}$must be long. Therefore we need only two shapes: ls⁺ and l⁺s to describethe intersection patterns of a line. ls⁺ applies when there are moreshort runs than long runs in the intersection pattern and l⁺s when thereare more long runs than short. In the case where there are equal numbersof runs, the two cases are equivalent.

The properties just described also apply to runs of runs. For example inintersection pattern 1004, the complete sequence of runs islslslsslslsslslss; therefore we have singularly occurring runs lseparating sequences of shorts runs s⁺ and the intersection pattern isconstructed of runs of runs with the shape ls⁺. The terminology of runsof runs is cumbersome so let us define these runs of runs to be secondorder runs. Where a run is defined by its position and length, a secondorder run is defined by its position, its length is defined by thenumber of runs which comprise it, and its shape is determined by whetherthere are more long or short runs in the pattern. By analogy with secondorder runs, we can define runs to be first order runs and pixels to bezero order runs.

To determine the possible lengths of the second order runs let usconsider the case where there are more short runs than long and continuethe use of our example. In this case, the number of second order runsmust be the same as the number of long runs in the pattern bydefinition. If there are n long runs in the line there are a runs todivide amongst n second order runs. Therefore, if the division is to beas even as possible the length of a short second order run, r^([2]),will be$r^{\lbrack 2\rbrack} = {\left\lfloor \frac{a}{n} \right\rfloor.}$There will be n^([2])≡a mod n long runs and a−n^([2]) short runs oforder 2. The second order runs appear at 1011 and 1013 in intersectionpattern 1004. There are more short runs 1013 than long, so the secondorder runs have the shape ls⁺, i.e., one long run followed by a sequenceof short runs. There are three long runs 1011 of length three and fourshort runs 1013 of length two.

There is no reason to end this hierarchical description at order 2. Wecan define a recursive hierarchical description of the intersectionpattern of the line l based on defining runs of higher order. For ordersthree and above, there is no restriction that the first run in theintersection pattern must be long and therefore all four possible shapescan occur. For order i, if there are more short runs of order i−1 in theintersection pattern, the shape of the order i runs will be s⁺l or ls⁺.If there are more long runs of order i−1, the shape will be l⁺s or ls⁺.An example of third order runs in the intersection pattern 1004 is shownat 1015 and 1017. There are three runs of order 3 amongst which theseven runs of order 2 are to be distributed as evenly as possible.Therefore the length of a short run of order 3 is${r^{\lbrack 3\rbrack} = {\left\lfloor \frac{7}{3} \right\rfloor = 2}},$and amongst the three order 3 runs, there will be 7 mod 3≡1 long run oflength r^([3])+1=3. The shape of the order 3 runs is s⁺l. For lines witha rational slope, the hierarchical description of ordered runs withinthe intersection pattern will be bounded, as the intersection pattern iseventually repeated if a and d are not coprime.

For the example line 1003 l: ${{l\text{:}y} = {\frac{17}{41}x}},$the process reaches its conclusion at order 4. There are three order 3runs, of which one is long and two are short. There will be thereforeone order 4 run containing all three order 3 runs and starting with theonly long order 3 run.

The order 4 run is the entire intersection pattern of the line l:${l\text{:}y} = {\frac{17}{41}x}$and the lattice

(41,17).Dealing with Lines that do not Intersect the Origin: FIG. 11

Introducing a non-zero intercept to a line with rational slope presentsa number of complications. Consider the line l:${l\text{:}y} = {{\frac{a}{d}x} + \beta}$where a and d are coprime. Within the frame of the lattice

(d,a), the line y forms an intersection pattern that is repeatedthroughout the infinite intersection pattern defined upon the unboundedlattice. Therefore we will only consider the intersection pattern

within the frame

(d,a).

Let us consider first the line l:${l\text{:}y} = {{\frac{a}{d}x} + \beta}$where β=0. Within the frame

(d,a), the vertical distance the line l is above any lattice pointwithin the intersection pattern is $\frac{b_{j}^{\lbrack 0\rbrack}}{d}$where b_(j) ^([0])=0, . . . , d−1 where the values b_(j) ^([0]) arethose values of the order 0 numerator sequence. The numerator sequenceis a sequence of the numerators of the fractions$\frac{b_{j}^{\lbrack 0\rbrack}}{d}.$These fractions specify the point in the left edge of the pixel at whichthe line intersects the pixel. Therefore, the line is at least$\frac{1}{d}$from any lattice point. If the value of β is raised slowly, theintersection pattern will remain unchanged until the line crosses orintersects a lattice point. Therefore, the first change will occur when$\beta = {\frac{1}{d}.}$Let the lattice point intersected be (x_(h), y_(h)).

As the line l: ${l^{\prime}\text{:}y} = {{\frac{a}{d}x} + \beta}$intersects the point (x_(h),y_(h)) and the line l:${l\text{:}y} = {\frac{a}{d}x}$intersects the origin, the intersection pattern

of the line l′ will be identical to the intersection pattern

of the line l translated by (x_(h), y_(h)). Therefore the introductionof an intercept to a line with rational slope will, more often than not,cause a shifting of the intersection pattern of the line. The key valuesof the intercept at which this translation of the intersection patternoccurs are $\beta = \frac{1}{d}$where b=0, . . . , d−1.

FIG. 11 gives an example. At 1001, the figure shows the intersectionpattern of the line l: ${l\text{:}y} = {\frac{17}{41}{x.}}$At 1101, the figure shows the effect of introducing an intercept valueof $\frac{23}{41}$on the interception pattern. The numerator sequence for order 0 of 1001is shown at 1109 and that for order 0 of 1101 at 1111.

Given that the numerator of the intercept is b=23, we have demarcatedthe pixels before (1103) and after (1104) b=23 in order 0 numeratorsequence 1109 by a dashed line 1102. The pixels 1104 to the right ofthis value are denoted by light gray and the pixels 1103 to the left bydark gray. At 1101, we can see that adding the intercept value$\frac{23}{41}$shifts the light gray pixels 1104 to the beginning of the intersectionpattern,

The dark gray pixels 1103 now form the end of the pattern

Coinciding with the shift in pixels, the values of numerator sequence oforder 0 1111 are also shifted and now start with a value of b₀^([0])=b=23. The shifting of the runs of all orders in the intersectionpattern with the introduction of a non-zero intercept is mimicked by theshift in the values of the numerator sequence.

The shifting of the intersection pattern due to the introduction of anon-zero intercept has a number of side effects. The initial and finalrun of any order may be truncated. This occurs when the numerator of theintercept is not in the numerator sequence of that order. A run order issplit and forms the initial and final partial run. If the numeratorsequence is to be calculated for this order, the initial numerator valuewill have to be calculated. For example, at 1101, the numerator value ofthe intercept is b=23. We know that all of the numerator values of order1 are less than p^([1])=a=17. Therefore the initial value of the order 0numerator sequence will not be the same as the initial value of theorder 1 numerator sequence and the initial and final runs of order 1will be truncated.

Using Hierarchies of Runs to Generate Lines

There are many line drawing techniques that take advantage of thestructure of a line's intersection pattern. At the pixel level, thestandard line drawing algorithm of Bresenham may be employed. In thisalgorithm, the point at which the line intersects the current pixeldetermines whether the next pixel's y coordinate is incremented by 1.See J. E. Bresenham, “An incremental algorithm for digital plotting”,ACM National Conference, August 1963. Bresenham's algorithm may be usedonly with lines whose start and end coordinates are rational numbers.Where the starting and end coordinates may be any real number, thewell-known DDA algorithm must be used. See A. Van Dam, J. Foley, S.Feiner and J. Hughes, Computer Graphics: Principles and Practice, SecondEdition in C, Addison-Wesley (1995). Like Bresenham's algorithm, DDAworks at the pixel level. Initially, floating-point x and y incrementsare computed. The x increment is the difference between the ending andstarting x coordinates divided by the line's length and the y incrementis the difference between the ending and starting y coordinates dividedby the line's length. Each time a pixel is set, the current x and ycoordinates, which are floating point values, are converted to integersand the pixel is set at the cell defined by the increment. Then the xincrement is added to the x coordinate and the y increment is added tothe y coordinate. A difficulty with any pixel-by-pixel approach is thatit requires a determination where the next pixel will be placed relativeto the last pixel for each new pixel. With the DDA algorithm, thisdetermination is a floating-point operation. As such, it is bothexpensive to perform and subject to rounding errors. Moreover, becausethe determination must be performed with every pixel, the roundingerrors may accumulate quickly and cause the line to be drawninaccurately.

The overhead of computing the location of the next pixel relative to thelast is avoided in algorithms that use runs of pixels instead ofindividual pixels to draw the line. At the level of a run of order 1,Reggiori has developed an algorithm that determines the length of thenext run in the line from the set of two possibilities. See G. B.Reggiori, “Digital computer transformations for irregular linedrawings”, Technical Report 403-22, New York University, April 1972.Stephenson generalizes these techniques to the full hierarchy of runs inthe line including runs of runs, runs of runs of runs, etc. See P.Stephenson, The structure of the digitized line, Ph.D thesis, James CookUniversity of North Queensland, 1998, which is incorporated by referenceherein. Like Bresenham's algorithm, the algorithms that use runs and runhierarchies are limited to lines whose start and end points haverational number coordinates.

All of these algorithms possess a similar conditional structureregardless of whether they are based on pixels such as Bresenham'spixel-based algorithm or the DDA algorithm, runs such as Reggiori'salgorithm, or a mixture of runs and runs of runs such as therun-length-slice algorithms. The slopes that are considered are boundedto lie in the range 0<α<1. For pixel-based algorithms, that limits thechoice of the next pixel to a possible set of two. For run-basedalgorithms, the choice is made between the two possible run lengths thatcan exist in the line. In all of these algorithms, the choice is made bychecking the value of a decision parameter against the value of zero.For values less than zero, one element of the possible set of choices isused; for values greater than zero, the other choice. For a value ofzero, each technique handles this case differently. For line drawingapplications, either choice is equally applicable and for ray tracing,this typically signifies a corner intersection.

Run-based line drawing may also be used to solve problems that arise inefficiently representing rays that traverse three-dimensional space andin speeding traversal of the three-dimensional space by such a ray. Itis an object of the present invention to provide solutions to suchproblems.

SUMMARY OF THE INVENTION

In one aspect, the invention is a technique employed in a computersystem for determining voxels in an object space that are intersected bya ray. Broadly stated, the technique involves making projections of theray on a plurality of planes in the object space, determining cells inthe planes that are intercected by the projections, and using theintersected cells to determine the intersected voxels. Determining whichcells are intercepted by the projections may be done using runs ofcells, as explained in the parents of the present application, and runsof order 1 or greater may be employed. The plurality of planes may be apair of intersecting planes which are perpendicular to each other andfor which the line of intersection is the ray's major axis.

Further details of this aspect include techniques for determiningintersected voxels that are 26-connected and intersected voxels that areonly six-connected. In the 26-connected case, the technique proceeds asfollows at the beginning of the next run of voxels:

-   -   if the first order runs of cells that include the points in the        projection corresponding to the beginning of the next run of        voxels end at the same coordinate of the ray's major axis, the        first order runs of cells together determine the next run of        voxels through the end of the first order runs; but    -   if the first order runs do not end at the same major axis        coordinate, the shorter first order run of cells and the        corresponding portion of the longer of the first order runs of        cells determine the next run of voxels.

In the 6-connected case, the above basic technique is used, but an extracell is added to the beginning of a first order run of cells prior tousing the first order run of cells to determine a run of voxels.

Further, this aspect includes techniques for determining from the runsof cells whether the voxels intersected by the ray have edge or cornerconnections. The voxels will have an edge connection if one of the firstorder runs of cells has a corner connection at a point and the otherdoes not; if both of the first order runs of cells have a cornerconnection at the point, the intersected voxels have a cornerconnection.

In another aspect, the invention is a technique employed in a computersystem of efficiently traversing a volume with a ray. The volume issubdivided into runs of voxels. Certain of the voxels in the runs areassociated with data that affects rays. Each ray that traverses thevolume intersects one or more of the first runs and it itself made up ofa set of second runs of voxels. The technique determines for a secondrun of voxels belonging to a particular ray whether the second runincludes a voxel of a first run that affects rays. When the second runincludes such a voxel, the data associated with the voxel is examined.It is advantageous for the technique if the volume has an axis which isthe major axis for both the particular ray and the first runs of voxelsand if there are three sets of first runs, each set having a differentaxis of the volume as its major axis. Continuing in more detail withthis aspect, the first runs contain significant runs that include voxelsthat affect rays. Determining whether the second run includes a voxelfrom a first run that affects rays includes determining whether thesecond run includes a voxel of a significant run.

In other details, aggregate information is associated with partitions ofthe first runs. The aggregate information indicates how one or more ofthe voxels in the partition affect rays. Determining whether the secondrun includes a voxel from a first run that affects ways includes usingthe aggregate information associated with a partition determine whetherthe partition contains a voxel that affects the particular ray. Thepartitions may be the significant runs. They may also be a number ofsets of partitions, with the partitions in each set having a differentlength. The partitions to use for a particular ray are selected inaccordance with the lengths of the runs in the ray.

Other objects and advantages will be apparent to those skilled in thearts to which the invention pertains upon perusal of the followingDetailed Description and drawing, wherein:

BRIEF DESCRIPTION OF THE DRAWING

FIG. 1 is an overview of the manner in which lines are drawn using thetechniques of the invention;

FIG. 2 is an overview of a system which draws lines using the techniquesof the invention;

FIG. 3 is a flowchart of how the bitmap processor processes completeruns of pixels;

FIG. 4 is a flowchart of how the bitmap processor processes truncatedruns of pixels;

FIG. 5 is a flowchart of how the bitmap processor processes a ran oforder 1;

FIG. 6 is a flowchart of how the first truncated run of an order i isprocessed and set;

FIG. 7 shows how the first pixel of a line is handled using thetechniques of the invention;

FIG. 8 shows a flowchart of how a complete run of order i is processed;

FIG. 9 shows a system in which the invention may be employed;

FIG. 10 shows the structure of a representation of a line as a set ofpixels;

FIG. 11 shows the effect of a displacement of the starting point fromthe origin on the structure of the representation;

FIG. 12 shows an example line;

FIG. 13 shows how the first pixel of the example line is handled;

FIG. 14 shows truncated runs of orders 1-3 in the example line;

FIG. 15 shows the geometry of the error term and structural parametersin a run of order 1;

FIG. 16 is a flowchart of a technique for drawing lines in parallel;

FIG. 17 is a diagram showing the geometry of ray traversal;

FIG. 18 is a diagram showing the geometry of the projection of a ray onplanes that are perpendicular to each other and intersect at the ray'smajor axis;

FIG. 19 is a first diagram showing how voxels intersected by a ray maybe computed from the cells intersected by the ray's projections;

FIG. 20 shows an algorithm for computing the voxels intersected by aray;

FIG. 21 is a second diagram showing how voxels intersected by a ray maybe computed;

FIG. 22 is a third such diagram showing how edge connections appear inthe projections;

FIG. 23 is an algorithm for detecting edge connections in the voxels;

FIG. 24 is a fourth diagram showing how corner connections appear in theprojections;

FIG. 25 is an algorithm for detecting corner connections in the voxels;

FIG. 26 is a diagram showing an encoding run in a volume;

FIG. 27 shows a ray traversing a volume which has been divided intoencoding runs;

FIG. 29 shows an algorithm for detecting an intersection between a rayand an encoding run;

FIG. 30 shows a ray traversing a volume in which the encoding runs havebeen divided into a hierarchy of runs; and

FIG. 31 shows an algorithm for using a hierarchy of runs to find voxelsthat will affect the ray.

DETAILED DESCRIPTION

As indicated in the section Cross references to related applications,this application contains the entire Detailed Description ofPCT/US02/24711 and the portion of PCT/US02/###### titled Paralleltechniques for determining the cells of a raster that are intercepted bya line. The material in the Detailed Description, which is new to thisapplication, begins with the section Problems of ray traversal in threedimensions.

Overview of Line Drawing with the Invention: FIG. 1

FIG. 1 is a flowchart [101] that provides an overview of how aline-segment may be drawn using this invention. Starting parametersinclude the maximum order of the runs to be used in drawing the lines[122] and the start and endpoints of the line [112,114]. The maximumorder can be predefined by the manufacturer of an apparatus based on theinvention, calculated based on the characteristics of the system orbased on the length of the line. In the preferred embodiment, themaximum order of the runs [120] is set to two.

The line segment defined between the starting point [112] (x₀, y₀) andthe end point [114] (x₁, y₁) consists of, at most, three sets of runs ofany order, i: The first truncated run of order i, the set of full-lengthruns of order i, and the final truncated run of order i. To draw theline segment from (x₀, y₀) to (x₁, y₁), we will describe how tocalculate and set in the Bitmap Memory [108]:

-   -   The first truncated run of order i [202].    -   The full length runs of order i [206].    -   The final truncated run of order i while setting the full length        runs of order i [206].

An alternative embodiment would calculate the first truncated run, eachfull length run and the final truncated run separately. Anotheralternative embodiment would have a desired maximum order [120] ofinfinity, and therefore the processing of the first truncated run lengthwould draw the entire line segment unless the maximum depth of thehierarchy of runs was reached.

Overview of the Components of a System for Drawing Lines According tothe Invention: FIG. 2

FIG. 2 shows the components of a system 201 for drawing lines accordingto the invention. The main components of system 201 are a Truncated RunProcessor [102], which calculates the structure of the first truncatedrun of order i [202] and the structure of the hierarchy of runs [140] inthe line segment, and the Run Processor [104], which calculates thelength of the full length runs and the final truncated run [206]. To setthe runs [202][208] into the Bitmap Memory [108], the description ofeach run is given to the Bitmap Processor [106]. The Bitmap Processor isresponsible for breaking the run into its composite pixels and settingeach pixel into the Bitmap Memory [108].

Truncated Run Processor [102]

The Truncated Run Processor [102] is responsible for calculating thestructure of the first truncated run of the desired order [120] andsupplying this structure to the Bitmap Processor [106] so that thetruncated run [194] can be set into the Bitmap Memory [108]. Bitmapmemory [108] may be a bit map such as bitmap 909 in FIG. 9.

While the Truncated Run Processor calculates the structure of the firsttruncated run, it also calculates the structure of the hierarchy of runs[140] and stores this information [142] for the use by the Run andMemory Processors [104,106]. During this process if the desired maximumorder of the process [120] is found to be greater than the maximum orderof the hierarchy of runs in the line, the maximum order [120] is set tobe the maximum depth of the hierarchy [124].

Based on the start and end points of the line [110], the Truncated RunProcessor also calculates and stores [132] the starting pixel positionfor the line [130].

Each component of the first truncated run is calculated based on theerror term [150]. As each component is processed, the error term isretrieved [152], updated and stored again [154]. Once the firsttruncated run has been processed, the error term remaining [150] is usedby the Run Processor [104] to define the full-length runs in a similarmanner.

Run Processor [104]

The Run Processor [104] is responsible for calculating the length ofeach full-length run of the desired maximum order in the digital line.These lengths are passed to the Bitmap Processor [106] such that theruns can be set into the Bitmap Memory [108]. The length of thefull-lengths is decided based on the error term [150], which isretrieved [156], updated and stored [158] to process the nextfull-length run. Thus, the error term is recalculated only at thebeginning of each full-length run of the maximum order instead of foreach pixel.

Bitmap Processor [106]

The Bitmap Processor [106] is responsible for actually setting the runsinto the Bitmap Memory [108].

We assume there is a method to set a pixel of the memory to a desiredvalue, similar to the setPixel (x, y, value) method described in theprior art. If a method to set a run of higher order, the BitmapProcessor [106] can take advantage of this function.

The Bitmap Processor [106] also keeps a track of the current position ofthe line being drawn in relation to the end point of the line [110]. Ifthe line moves reaches the end point [114], the Bitmap Processor [106]signals the Truncated Run Processor [102] or the Run Processor [104] toterminate [162] using an internally stored termination condition [160].

Bitmap Memory [108]

The Bitmap Memory [108] is the memory representing the raster display ofa graphics system, video memory, image memory or any other type ofbitmap.

Details of Operation of the Components of System 201: FIGS. 3-8

For the remainder of this description we will describe:

-   -   1. How a truncated or full-length run of order i [202][208] is        set into the Bitmap Memory [108] via the Bitmap Processor [106].    -   2. How the first truncated run of order i in the line segment is        processed [202] by the Truncated Run Processor [102].    -   3. How the full-length runs of order i in the line segment [208]        are processed by the Run Processor [104].        How a Run of Order i is set into the Bitmap Memory via the        Bitmap Processor: FIG. 3

The Bitmap Processor [106] handles the process of setting the pixels ofthe Bitmap Memory [108]. To draw the line, the Bitmap Processor is givena set of commands [194][196] to draw a collection of runs of variousorders into the Bitmap Memory. The information given in the command todraw a run is the length of the run, r_(j) ^([i]), and the order of therun, i [194][196]. The position of the run in the bitmap [130], (x_(j)^([i]), y_(j) ^([i])), is initialized [132] by the Truncated RunProcessor [102] and is retrieved by the Bitmap Processor [136]. As therun is drawn, the coordinate is updated to the position of the next runof order i, and is stored for the next iteration [134].

The process of setting a run in the Bitmap Memory, involves reducing therun to a set of pixels and setting each individual pixel in the BitmapMemory. Reducing the run of order i into runs of order 0 (pixels) isperformed by recursively reducing the run of order k to runs of orderk−1 for each order from k=i to k=2 resulting in a description of the runof order i in terms of runs of order 1. Each pixel in the run of order 1is then set into memory.

We therefore describe:

-   -   How a full-length run of order i is reduced into runs of order        i−1.    -   How the first truncated run of order i is reduced into runs of        order i−1.    -   How the pixels of a run of order 1 are set into the Bitmap        Memory.

In an alternative embodiment, if the Bitmap Memory [108] permits thesetting of runs of pixels or runs of runs, etc., the process of reducingthe run of order i to runs of order i−1 can be stopped at the permittedlevel. Therefore, instead of setting pixels into the Bitmap Memory, runsor runs of runs can also be used.

Reducing a Full-length Run of order k into Runs of Order k−1: FIG. 3

As described in the Description of related art, each run of any nonzeroorder i within the digital line has, at most, only two possible runlengths that are consecutive integers where the length is measured asthe number of composite runs of the next lowest order, i−1. As aconsequence of this, we can denote the occurrence of a short run by thesymbol s and a long run by the symbol l. For example the runs of order 1in the line segment (0, 0) to (41, 17) can be described by(lslslsslslsslslss)^([1]).

A run of order i can then be defined recursively by:

-   -   A run of order 1 is the maximal set of contiguous pixels with a        similar abscissa.    -   A run of order i is the maximal set of contiguous runs of order        i−1 such that each run has a shape (l⁺s)^([i−1]), (ls⁺)^([i−1]),        (s⁺l)^([i−1l) , or (sl⁺)^([i−1]), where (l⁺)^([i−1]) and        (s⁺)^([i−1]) denote the occurrence of one or more long and short        runs of order i−1 respectively.

We note that a run of order 2 can only possess the shape (l⁺s)^([1]) or(ls⁺)^([1]).

Each run of any given order in the line segment has the same shape. Todetermine for a given order which of the four shapes is to be used, wedefine the runs of each order to have a type, t^([i]), which has a value0 or 1. For a run of order i:

-   -   If t^([i−1])=0 and t^([i])=0 the shape of the run of order i is        s^([i])=0: l⁺s.    -   If t^([i−1])=0 and t^([i])=1 the shape of the run of order i is        s^([i])=1: ls⁺.    -   If t^([i−1])=1 and t^([i])=0 the shape of the run of order i is        s^([i])=2: sl⁺.    -   If t^([i−1])=1 and t^([i])=1 the shape of the run of order i is        s^([i])=3: s⁺l.

The type of order 1, t^([1]), is defined to be zero. The types of ordersgreater than 1 are calculated using structural parameters; thestructural parameters and their calculation will be described in detailbelow.

To reduce a run of order i to the composite runs of order i−1, theprocedure 301 shown in FIG. 3 is followed.

Inputs to the process are the length of the run of order i to be set[194] and the structure of the hierarchy of runs. The process will ceaseif the end point of the line has been reached [304].

If the order of the run is 1 [306], the pixels of the run are setdirectly [340] into the Bitmap Memory.

If the order of the run is at least 2, the run is reduced into itscomposite runs of the next lesser order and each run is set into theBitmap Memory:

-   -   If t^([i−1])=0 and t^([i])=0 the shape of the run is s^([i])=0:        l⁺s, therefore we set in order r_(j) ^([i])−1 long runs of order        i−1 (length r_(l) ^(i−1])=r_(s) ^([i−1)+1) and one short run of        order i−1 (length r_(s) ^([i−1])).

This process follows the path from [308], [310], [312], [314].

If t^([i−1])=0 and t^([i])=1 the shape of the run is s^([i])=1: ls⁺, weset one long run of order i−1 (length r_(l) ^([i−1])=r_(s) ^([i−1])+1)and r_(j) ^([i])−1 short runs of order i−i (length r_(s) ^([i−1])).

This process follows the path from [308], [310], [316], [318].

-   -   If t^([i−1])=1 and t^([i])=0 the shape of the run is s^([i])=2:        sl⁺, we set one short run of order i−1 (length r_(s) ^([i−1]))        and r_(j) ^([i])−1 long runs of order i−1 (length r_(l)        ^([i−1])=+1).

This process follows the path from [308], [320], [322], [324].

-   -   If t^(i−1])=1 and t^([i])=1 the shape of the run is S^([i])=3:        s⁺l, we set r_(j) ^([1])−1 short runs of order i−1 (length r_(s)        ^([i−1])) and one long run of order i−1 (length r_(l)        ^([i−1])=r_(s) ^([i−1]) +1).

This process follows the path from [308], [320], [326], [328].

The type and length of the short and long runs of any order are storedin the structure of the hierarchy of runs [146].

Reducing the First Truncated Run of Order k into Runs of Order k−1: FIG.4

The process of reducing the first truncated run of order k intofull-length runs of order k−1 follows directly from the instructions toreduce full-length runs. The only variation comes from the fact that atruncated run is defined to be a run not of its full-length. Thereforeat least one run of order k−1 has been truncated from the run of orderk. Because this is the first run in the line, the truncated run occursat the beginning of the run. In a preferred embodiment, the truncatedruns received by the Bitmap Processor have a length of at least one.Therefore the only situations that need to be handled differently arewhen shapes s^([i])=1: ls⁺ and s^([i])=2: sl⁺ occur.

The process is shown at 401 in FIG. 4. As with the process for afull-length run, the inputs are the length of the run of order i that isto be set [194] and the structure of the hierarchy of runs [146]. Theprocess 401 will cease if the end point of the line has been reached[404].

If the order of the run is 1 [406], the pixels of the run are setdirectly [440] into the Bitmap Memory.

If the order of the run is at least 2, the run is reduced into itscomposite runs of the next lesser order and each run is set into theBitmap Memory:

-   -   If t^([i−1])=0 [408] and t^([i])=0 [410], the shape of the run        is s^([i])=0: l⁺s. As the run is truncated, at least one run of        order i−1 is removed from the beginning of the run. As the        length of the truncated run is r₀ ^([i]), this leaves r₀        ^([i])−1 long runs of order i−1 and the final short run of order        i−1. For example if the first run at full length is the sequence        of runs (lllls)^([i−1]) and the initial truncated run length is        r₀ ^([i])=3, the truncated run comprises two long runs and the        final short run of order i−1: (lls)^([i−1]). Therefore to set        the initial truncated run length, we set in order r₀ ^([i])−1        long runs of order i−1 (length r_(l) ^([i−1])=r_(s) ^([i−1])+1)        [412] and one short run of order i−1 (length r_(s) ^([i−1]))        [414].    -   Similarly, if t^([i−1])=0 [408] and t^([i])=1 [410], the shape        of the run is s^([i])=1: ls⁺. The truncated run has less runs        than a short run of order i, therefore given the shape of the        run ls+, there can only be short runs of order i−1 comprises the        truncated run. Therefore we set r₀ ^([i]) short runs of order        i−1 (length r_(s) ^([i−1])) [418].    -   If t^([i−1])=1 [408] and t^([i])=0 [420] the shape of the run is        s^([i])=2: sl⁺, we set r₀ ^([i]) long runs of order i−1 (length        r_(l) ^([i−1])=r_(s) ^([i−1])+1) [424].    -   If t^([i−1])=1 [408] and t^([i])=1 [420] the shape of the run is        s^([i])3: s⁺l, we set r₀ ^([i])−1 short runs of order i−1        (length r_(S) ^([i−1])) [426] and one long run of order i−1        (length r_(l) ^([i−1])=r_(s) ^([i−1])+1) [428].

The type and length of the short and long runs of any order are storedin the structure of the hierarchy of runs [146].

Setting the Pixels of a Run of Order 1 into the Bitmap Memory: FIG. 5

FIG. 5 is a flowchart 501 showing how to set the pixels of a run oforder 1 [340]. The method requires as inputs the starting pixel [136],(x_(j) ^([1]), y_(j) ^([1])), the end point of the line [114], and thelength of the run [194], r_(j) ^([1])), are required. The starting pixelis stored and updated internally [136]. The length of the run [194] ispassed to the Bitmap Processor.

Each pixel in the run has a similar y-coordinate, which is they-coordinate of the starting pixel [136], y_(j) ^([1]).

If the y-coordinate for the run being set is the same as the floor ofthe y-coordinate of the end point of the line, y_(j) ^([1])=└y₁┘, therun being set is the last run in the line [402]. In this case:

-   -   The termination condition [160] signaling that the end point of        the line has been reached is set [406].    -   The length of the run to be set [194] is changed to the        difference between the x-coordinate of the start of the run        being set and the ceiling of the x-coordinate of the endpoint of        the line, r_(j) ^([1])=┌x₁┐−x_(j) ^([1]) [404].    -   The run of pixels of length r_(j) ^([1])=┌x₁┐−x_(j) ^([1]) is        set in the Bitmap Memory at location (x_(j) ^([1]), y_(j)        ^([1])).

To set each of the pixels in the run into the Bitmap Memory, for each ofthe pixels in the run [410][420][426] starting at the pixel (x_(k)^([0]), y_(k) ^([0]))=(x_(j) ^([1]), y_(j) ^([1])) [412]:

-   -   Set the pixel (x_(k) ^([0]), y_(k) ^([0])) into the Bitmap        Memory [422].    -   Move to the position of the next pixel, $\begin{matrix}        {\left( {x_{k + 1}^{\lbrack 0\rbrack},y_{k + 1}^{\lbrack 0\rbrack}} \right) = {\left( {{x_{k}^{\lbrack 0\rbrack} + 1},y_{k}^{\lbrack 0\rbrack}} \right).}} & \lbrack 424\rbrack        \end{matrix}$

Once every pixel in the run has been set into the Bitmap Memory, thestarting coordinate of the next run of order 1 in the line is calculatedfrom the current coordinate value by incrementing the currenty-coordinate of the next pixel position(x_(j + 1)^([1]), y_(j + 1)^([1])) = (x_(k)^([0]), y_(k)^([0]) + 1).How the First Truncated Run of Order i in the Line Segment is Processedby the Truncated Run Processor: FIGS. 6 and 12

To calculate the length of the initial truncated run length of themaximum order i [120] we calculate and set in the Bitmap Memory eachinitial truncated run length of order 1 to i in turn. This is shown inflowchart 601 of FIG. 6. If no truncated run length exists for anyorder, none is set. The process stops if:

-   -   The endpoint of the line [110] is reached [514].    -   The maximum depth of the hierarchy is reached [520].    -   The desired maximum order is reached [530].

A truncated run length of order k is of course made up of run lengths oforders 0 through k−1. In the above technique, the orders 1 through k areprocessed beginning with order i=1 and moving order by order to i=k.With each of the orders, the bits for that order are set before the nexthigher order is processed. For a given order i, the bits set for orderi−1 will always make up the complete truncated or untruncated first runof the order i.

Initialize the Slope of the Line [502].

The endpoints of the line segment are the real number coordinates (x₀,y₀) and (x₁, y₁) [112], such that x₀≦x₁ and y₀≦Y₁, therefore the slopeof the line, α, has a value between zero and one, 0≦α<1.

The slope of the line is calculated by α=(y₁−y₀)/(x₁−x₀) [502]. Inexample line 1201 of FIG. 12, the equation of the line is y=17x/41+7/82,therefore α=17/41. Example line 1201's endpoint (x₀, y₀) 1203 is(1/2,12/41); its endpoint (x₁, y₁) 1207 is (41+1/2,17+12/41). Note thatwhile the coordinates of the endpoints are rational numbers in thiscase, a coordinate may be any real number.

Initialize Order 1 Starting Pixel and Error Term [504]: FIGS. 7 and 13

The first step in the algorithm is to handle the pixel in which thestart point of line lies as described in FIG. 13. At the completion ofthis step, the dark pixel 1301 at the start of the line will have beenprocessed.

As may be seen from FIG. 13, the coordinates of the lower left-handcorner 1303 of pixel 1301 that contains the start point 1203 of line1201 is:x ₀ ^([0])=└1/2┘=0y ₀ ^([0])=└12/41┘=0

To decide whether to set the first pixel in the line separately from thefirst run we define the distance X₀ ^([0]) 703 from the start point ofthe line to the beginning of the pixel as described in FIG. 7:X ₀ ^([0]) =x ₀ −x ₀ ^([0])or 1/2 in our example. If x₀ ^([0])≠0, the first pixel must be setseparately from the first run.

We also require the initial value error term {circumflex over (β)}₀^([0]) of order zero normalized against the slope of the line. Given theaverage run length in the line τ^([1])=41/17 in example line 1201,{circumflex over (β)}₀ ^([0])=(y ₀ =y ₀ ^([0]))τ^([1])=12/17where y₀−y₀ ^([0]) is the distance from the y coordinate of startingpoint 203 to the y coordinate of the lower left-hand corner 1303 ofpixel 1301 and $\tau^{\lbrack 1\rbrack} = {\frac{1}{\alpha^{1}}.}$Since X₀ ^([0]) is nonzero, the Bitmap Processor sets a run of order 0(a pixel) into the bitmap memory at the location (x₀ ^([0]), y₀^([0]))=(0,0). There are two choices for the position of the next pixelin the line segment (1,0) or (1,1). The choice from these two is madeusing the value of the error term {circumflex over (β)}₀ ^([0]). Thereare two choices for the position of the next pixel in the line segment(x₀ ^([0])+1, y₀ ^([0])) and (x₀ ^([0])+1, y₀ ^([0])+1). The choice fromthese two is made by the value of the error term {circumflex over (β)}₀^([0]). Firstly, increment the x coordinate and calculate the next errorterm of order 0:x ₀ ^([0]) =x ₀ ^([0])+1{circumflex over (β)}₀ ^([0])={circumflex over (β)}₀ ^([0])+1−X ₀ ^([0])In the example of FIG. 7, this comes out to 12/17+1−1/2=41/34, which isless than τ¹, which equals 41/17.

If the error term is less than the average run length of order 1 asshown at 701 in FIG. 7, the next pixel is pixel (x₀ ^([0])+1, y₀^([0])).if {circumflex over (β)}₀ ^([0])<τ^([1]) theny₀ ^([1])=y₀ ^([0])

If the error term is greater than or equal to the average run length oforder 1 as shown at 709 in FIG. 7, the next pixel is pixel (x₀ ^([0])+1,y₀ ^([0])+1). Accordingly, the error term is also updated.if {circumflex over (β)}₀ ^([0])≧τ^([1]) theny₀ ^([0]) =y ₀ ^([0])+1{circumflex over (β)}₀ ^([0])={circumflex over (β)}₀ ^([0])−τ^([1])

If X₀ ^([0]) is zero, the values of position of the pixel and the errorterm are not altered and no pixel is set into the memory.

The position of the next run of order 1 in the line segment and thefirst error term of order 1 is now the position of the next order 0 runand error term of order 0.x₀ ^([1])=x₀ ^([0])y₀ ^([1])=y₀ ^([0])

The coordinate (x₀ ^([1]), y₀ ^([1])) is stored as the position of thenext run in the line segment [130].

The error term {circumflex over (β)}₀ ^([0]) is stored [152] as theerror term for the next run [150].

Setting the First Truncated Run of Order i>0

Once the first order has been initialized [504][506], the firsttruncated run of order 1 through i is set, where i is either the desiredmaximum order [122] or the maximum order in the hierarchy of runs in theline [124].

The process for setting the first truncated run of order k=1 is asfollows. Starting at order k=1 [506]:

-   1. Calculate the length of the first truncated run of order k [510],    r₀ ^([k]).-   2. Set the truncated run of order k into the Bitmap Memory [512] if    the run is truncated.-   3. Check that the end point of the line has not been reached [514].    If it has, cease processing and if it has not continue.-   4. Calculate the next error term of order k [516].-   5. Check if k is the maximum depth of the hierarchy of runs in the    line [520]. If so, set the maximum order of runs [120] to be the    maximum order of the hierarchy of runs [124], update and store [152]    the final value of the error term [150] [590] and cease processing.    If not continue.-   6. Check if k is the desired maximum order of runs for the line    [530]. If so, update and store [152] the final value of the error    term [150] [590] and cease processing. If not continue.-   7. The current error term of order k is the next error term for    order k+1 [532].-   8. Move to the next order [534], k=k+1.    Structure of the Hierarchy of Runs

During this process the structure of the hierarchy of runs will becalculated and stored [140]. The structure of the hierarchy of runsstores a description of each level in the hierarchy from order 1 to thedefined maximum order [120]. Each record of the description for order istores the following entries:The slope of the line of order i, α^([k]).

-   -   -   The slope of the line of order 1 is defined to be the slope            of the line, α^([1])=α.        -   The slope of the line of order k>1, αk^([k])=min({circumflex            over (μ)}^([k−]), {circumflex over (v)}^([k−1])).            The average length of a run of order i in the line,            τ^([k])=1/α^([k]).            The length of a short run of order i, r_(s)            ^([k])=└τ^([k])┘.            The length of a long run of order i, r_(l) ^([k])=|τ^([k])|.            The type of the runs of order i in the line, t^([k]), either            0 or 1.        -   The type of order 1 is defined to be t^([1])=0.        -   The type of order k>1 is defined to be t^([k])=0 if            {circumflex over (μ)}^([k−1])≦{circumflex over (v)}^([k−1]).        -   The type of order k>1 is defined to be t^([k])=1 if            {circumflex over (μ)}^([k−1])>{circumflex over            (v)}^([k−1])).            The structural parameters {circumflex over (μ)}^([i]) and            {circumflex over (v)}^([i]).        -   If the type of order k, t^([k])=0, then the structural            parameters of order k,            {circumflex over (v)}^([k])=τ^([k]) −r _(s) ^([k]), and            {circumflex over (μ)}^([k])=1−{circumflex over (v)} ^([k])        -   If the type of order k, t^([k])=1, then the structural            parameters of order k,            {circumflex over (μ)}^([k])=τ^([k]) −r _(s) ^([k]), and            {circumflex over (v)}^([k])=1−{circumflex over (μ)}^([k])            Geometry of the Error Parameter {circumflex over (β)}_(j)            ^([1]) and the Structural Parameters {circumflex over            (υ)}^([1]) and {circumflex over (μ)}^([1]): FIG. 15

As may be seen from the foregoing, the error parameter {circumflex over(β)}_(j) ^([i]) is used to calculate lengths of runs and the structuralparameters {circumflex over (υ)}^([i]) and {circumflex over (μ)}^([i])are used to determine the type of an order. FIG. 15 shows geometricallywhy the error parameter and the structural parameters can be used inthis fashion. The figure shows a first-order run j of pixels 1505representing a portion of line 1507 with slope α. In first-order runs,the slope αl^([1]) of the run is the same as the slope of the line 1507.Shown at 1503 is the last pixel of the run j−1. The coordinates (x_(j)^([1]), y_(j) ^([1])) are the coordinates of the lower left-hand cornerof the first pixel of run j; the coordinates(x_(j + 1)^([1]), y_(j + 1)^([1]))are the coordinates of the lower left-hand corner of the first pixel ofthe first pixel of run j+1. The values of the error parameter and thestructural parameters are defined geometrically by means of lines thatare parallel to line 1507 and intersect corners of the pixels of run j.Thus, line 1509 intersects the upper left-hand corner of pixel 1503;line 1511 intersects the upper left-hand corner of the last pixel 1506in run j; line 1513 intersects the bottom left-hand corner of the firstpixel in run j.

The error parameter {circumflex over (β)}_(j) ^([1]) is the distance onthe line formed by the top of pixel 1506 between the intersections oflines 1507 and 1513 with that line. The greater the distance betweenlines 1507 and 1513, the shorter the next run must be. Moreover, since aline is drawn beginning with a partial run of the maximum order k andthe partial run is made by beginning with a partial run of order 1,followed by a partial run of order 2, and continuing through order k,the error term for the partial run of order i is that resulting fromorder i−1.

The structural parameter {circumflex over (υ)}^([1]) is the distance onthe line formed by the top of pixel 1506 between the intersection ofline 1513 with the line formed by the top of pixel 1506 and the upperleft-hand corner of pixel 1506. The type of order 1 is defined to be 1,so the structural parameter {circumflex over (μ)}^([1])=1−{circumflexover (υ)}^([1]). The structural parameters vary with the slope α^([i])of the line as represented by runs of order i, and thus can be used todetermine the type of order i+1.

Calculation of the Length of the Truncated Run

To calculate the length of the first truncated run of order k [510], r₀^([k]), first check if the run is truncated. If the run is nottruncated, the length of the first run is assumed to be zero, r₀^([k])=0, and no run is set into the Bitmap Memory [512]. The first runof order k in the line is truncated if the stored error term [150],{circumflex over (β)}₀ ^([k−1])>1.

If the run of order k is truncated, the length of the run is:r₀ ^([k])=|τ^([k−1])| if the type t^([k])=0.r₀ ^([k])=|{circumflex over (β)}₀ ^([k−1])| if the type t ^([k])=1.

If the run is not truncated, the error term for the next order is thesame as this order {circumflex over (β)}₀ ^([k])={circumflex over (β)}₀^([k−1]). The run of order k and length r₀ ^([k]) can now be set intothe Bitmap memory.

If the run is truncated the initial truncated run length and initialdecision value of order k is:{circumflex over (β)}₀ ^([k]) =r ₀ ^([k])−τ^([k])+{circumflex over(β)}^([k−1]) if the type t ^([k])=0.{circumflex over (β)}₀ ^([k−1]) −r ₀ ^([k]) if the type t ^([k])=1.

At the end of the process of setting the initial truncated runs, we havedrawn the first truncated run of order i in the line segment. We havethe position of the next run of order i in the line if it exists, (x₁^([i]), y₁ ^([i])) and the decision value {circumflex over (β)}₁ ^([1]).We therefore move to the next phase of the process which calculates thelength of the next run of order i in the line, which is then set intothe memory. The length of the next run of order i is decided from thetwo possible by the value of the decision variable, {circumflex over(β)}₁ ^([i]).

To update [590] the value of the error term of order i, {circumflex over(β)}₁ ^([i]), the value of {circumflex over (β)}₁ ^([i]) is adjusted by−{circumflex over (v)}^([i]), {circumflex over (β)}₁ ^([i])={circumflexover (β)}₁ ^([i])−{circumflex over (v)}^([i]) and this value is stored[152] for use by the Run Processor [104].

An Example of Setting Orders 1-3: FIG. 14

FIG. 14 shows line 1201 and the pixels 1402 that will be generated torepresent it according to the technique under discussion. At 1401 isshown how the first truncated run 1403 of order 1 is set; at 1405 isshown how the first truncated run 1405 of order 2 is set; and at 1407how the first truncated run of order 3 is set.

Once the first order has been initialized [504][506], the firsttruncated run of order 1 through i is set, where i is either the desiredmaximum order [122] or the maximum order in the hierarchy of runs in theline [124]. For the sake of the example, we will take the desired orderto be i=3. As part of the initialization, first pixel 1301 has been setand the position of the next pixel relative to first pixel 1301 has beendetermined as described above.

Order 1

The next step is to handle the first truncated run of order 1 if itexists. In the example line, the portion on the line that will beprocessed is described by the dark pixels 1403.

The definition of the hierarchy of runs at order 1 in [146] is:The slope of the line of order 1, α^([1])=17/41.The average length of a run of order 1 in the line, τ^([1])=41/17.The length of a short run of order 1, r_(s) ^([1]=└τ) ^([1])=2.The length of a long run of order 1, r_(l) ^([1])=|τ^([1])|=3.The type of the runs of order 1 in the line is defined to be, t^([1])=0.The structural parameters {circumflex over (μ)}^([1]) and {circumflexover (v)}^([1]).{circumflex over (v)}^([1])=τ^([1]) −r _(s) ^([1])=7/17, and{circumflex over (μ)}^([1])=1−{circumflex over (v)} ^([1])=10/17.

The length of the first truncated run of order 1 [510], r₀ ^([1]), istruncated since the stored error term [150], {circumflex over (β)}₀^([0])=41/34>1. The length of the truncated run is r₀^([1])=|τ^([1])−{circumflex over (β)}₀ ^([0])|=┌41/17−41/34┐=2 since thetype t^([1])=0. The run of order 1 and length r₀ ^([1])=2 can now be setinto the Bitmap memory at position (x₀ ^([0]), y₀ ^([0]))=(1,0). Theinitial decision value of order 1 is {circumflex over (β)}₀ ^([1])=r₀^([1]−τ) ^([1])+{circumflex over (β)}₀ ^([0])=2−41/17+41/34=27/34.

Order 2

The initial truncated run of order 2 1405 is described by the darkpixels 1407 in line 1201. The initial truncated run of order 2 has alength of one (i.e. one run of order 1). The position of the initialtruncated run of order 2 is (x₀ ^([2]), y₀ ^([2]))=(x₀ ^([1]) +r ₀^([1]), y₀ ^([2])+1)=(3,1).

The hierarchy of runs at order 2 is described by the parameters:

-   -   a The slope of order 2, α^([2])=min({circumflex over (μ)}^([1]),        {circumflex over (v)}^([1]))=7/17.    -   The average length of a run of order 2 in the line,        τ^([2])=17/7.    -   The length of a short run of order 2, r_(s) ^([2])=2.    -   The length of a long run of order 2, r_(l) ^([2])=3.        As {circumflex over (μ)}^([1])>{circumflex over (v)}^([1]):        t^([2])=1        {circumflex over (μ)}^([2])=τ^([2]) −r _(s) ^([2])=3/7        {circumflex over (v)}^([2])=−{circumflex over (μ)}^([2])=4/7

The first run of order 2 in the line is truncated because the storederror term [150] renormalized to the slope of order 2, {circumflex over(β)}^([1])={circumflex over (β)}₀ ^([1])τ^([2])=27/14, is greaterthan 1. Therefore as t^([2])=1:r ₀ ^([2])=|{circumflex over (β)}₀ ^([1])|=1{circumflex over (β)}₀ ^([2])={circumflex over (β)}₀ ^([1]) −r ₀^([2])=13/14

As the shape of the runs of order 2 is shape ls⁺, the last run of order1 in a run of order 2 is a short run of order 1. A short run of order 1has a length of two pixels. Therefore two pixels are set in this stage.

Order 3

The initial truncated run of order 3 1409 is described by the darkpixels 1411. The initial truncated run of order 3 has a length of one(i.e. one run of order 2). The runs of order 3 in the line have a shapes⁺l, the last run of order 2 in a run of order 3 is a long run.Therefore we will set a run of order 2 of length three.

The hierarchy of runs at order 3 is described by the parameters:

-   -   The slope of order 3, α^([3])=min({circumflex over (μ)}^([2]),        {circumflex over (v)}^([2]))=3/7.    -   The average length of a run of order 3 in the line, τ^([3])=7/3.    -   The length of a short run of order 3, r_(s) ^([3])=2.    -   The length of a long run of order 3, r_(l) ^([3])=3.        As {circumflex over (μ)}^([2])<{circumflex over (v)}^([2]):        t^([3])=0        {circumflex over (v)}^([3])=τ^([3]) −r _(s) ^([3])=2/3        {circumflex over (μ)}^([3])=1−{circumflex over (v)} ^([3])=2/3

The first run of order 2 in the line is truncated because the storederror term [150] renormalized to the slope of order 3, {circumflex over(β)}₀ ^([2])={circumflex over (β)}₀ ^([2])τ^([3])=13/6, is greaterthan 1. Therefore as t^([2])=0:r₀ ^([3])=|τ^([3])−{circumflex over (β)}₀ ^([2])|=┌1/6┐=1{circumflex over (β)}₀ ^([3]) =r ₀ ^([3])−τ^([2])={circumflex over (β)}₀^([2])=5/6

Therefore, a single run of order 2 comprises the run of order 3. Nowthat the order chosen for the iterative portion of the technique hasbeen reached, the normalized intercept value of order 3, {circumflexover (β)}₀ ^([3])=5/6, can be used to initialize the iterative decisionprocess.

Setting Full-length Runs of Order i in the Line; FIG. 8

To set a full length run of order i into the Bitmap Memory [108], theRun Processor [104] performs the steps shown in flowchart 801. Theinputs are the maximum order of the runs (120), the structure of thehierarchy of runs [146] which was determined by truncated run processor[102] while processing the truncated runs, and the error term{circumflex over (β)}_(j) ^([i]) provided by truncated run processor[102]. For each run of order i, the run processor determines the lengthof the run of order i in the line and passes the run's length and orderto the Bitmap Processor [106] to set the actual pixels of the run intothe Bitmap Memory [732]. The calculation of the length of the run oforder i is based on the type of order (i) [702] and the value of theerror term [150], {circumflex over (β)}_(j) ^([i]), retrieved [158] bythe Run Processor [710,720]. Once the length of the run has beencalculated, the error term is updated, β̂_(j + 1)^([i]),and stored [156]. The flow chart branches involved here are [702, 710,712, 714]; [702, 710, 716, 718]; [702, 720, 722, 724]; and [702, 720,726, 728].

The calculation of the length of the run of order i in the line, r_(j)^([i]) is as follows:If type t^([i])=0:

-   -   -   If {circumflex over (β)}_(j) ^([i])≧0 [710] the run is short            [712], r_(j) ^([i])=r_(s) ^([i]), and $\begin{matrix}            {{\hat{\beta}}_{j + 1}^{\lbrack i\rbrack} = {{\hat{\beta}}_{j}^{\lbrack i\rbrack} - {{\hat{v}}^{\lbrack i\rbrack}.}}} & \lbrack 714\rbrack            \end{matrix}$        -   If {circumflex over (β)}_(j) ^([i])<0 [710] the run is long            [716], r_(j) ^([i])=r_(l) ^([i]), and $\begin{matrix}            {{\hat{\beta}}_{j + 1}^{\lbrack i\rbrack} = {{\hat{\beta}}_{j}^{\lbrack i\rbrack} + {{\hat{\mu}}^{\lbrack i\rbrack}.}}} & \lbrack 718\rbrack            \end{matrix}$            If type t^([i])=1:        -   If {circumflex over (β)}_(j) ^([i])≧0 [720] the run is long            [722], r_(j) ^([i])=r_(j) ^([i]), and $\begin{matrix}            {{\hat{\beta}}_{j + 1}^{\lbrack i\rbrack} = {{\hat{\beta}}_{j}^{\lbrack i\rbrack} - {{\hat{v}}^{\lbrack i\rbrack}.}}} & \lbrack 724\rbrack            \end{matrix}$        -   If {circumflex over (β)}_(j) ^([i])<0 [720] the run is short            [726], r_(j) ^([i])=r_(s) ^([i]), and $\begin{matrix}            {{\hat{\beta}}_{j + 1}^{\lbrack i\rbrack} = {{\hat{\beta}}_{j}^{\lbrack i\rbrack} + {{\hat{\mu}}^{\lbrack i\rbrack}.}}} & \lbrack 728\rbrack            \end{matrix}$

The values of the structural parameters {circumflex over (μ)}^([i]) and{circumflex over (v)}^([i]) and the length of a short and long run,r_(s) ^([i]) and r_(l) ^([i]), are retrieved [144] from the structure ofthe hierarchy of runs [140].

Once the length of the run is calculated, the run is set into memory[194] by the Bitmap Processor [106] and [732]. The Run Processor keepsprocessing runs as just described until it reaches the end of the linesegment [730].

Applications of the Line Drawing Techniques

The line drawing techniques described herein can be employed in anyapplication where techniques of the general class to which thesetechniques belong are useful. The techniques are of course particularlyuseful in applications where the coordinates of the endpoints of thelines may have real number values, including values that are irrationalnumbers. Such applications include the following:

Two-Dimensional Polygon Filling and Scan Conversion

Techniques for filling a polygon on a raster or bitmap display oftenstep along the circumference of the polygon and fill or scan convert thepolygon along the horizontal runs defining the interior of the polygon(see Lathrop 1990). A problem with using existing line drawingtechniques for this purpose is that the endpoints of the line aredefined at pixel or sub-pixels positions and using these algorithms tostep along the circumference of the polygon can cause dropouts ormultiply processed pixels to occur. The line drawing technique describedcan be used to step along the circumference of the polygon without theseerrors.

Three-Dimensional Polygon Scan Conversion

When a 3D polygon is being rendered into a 2D image using interpolativetechniques such as Gouraud or Phong shading, typically a techniquesimilar to the 2D polygon scan conversion embodiment is used. Atechnique that is also possible is to shade the polygon using lines ofconstant depth. Think of a plane parallel to the view plane being used.As the plane is moved backwards through the polygon, the line ofintersection between the plane and the polygon can be defined. This linehas a constant depth that can be used to make the shading technique moreefficient. To shade the polygon you can step along this line using theline drawing technique described setting each pixel the color calculatedby the shading algorithm.

Visibility Checking—Computer Games and Simulation In computer games andsimulation a common problem is to determine if two points are visible toone another. If the environment that can occlude the two points can bedefined on a 2D grid, the line drawing technique describe can be used tocheck each of the grid points linearly separating the two points todetermine if they are occluded from one another.

Visibility Checking—Mobile Network Planning

In the mobile phone industry, one of the greatest problems is inplanning the mobile phone base station network to ensure effectivecoverage. Mistakes in deciding where to place base stations and theirpower requirements can be very expensive. One common solution is to usetwo-dimensional ray tracing methods to simulate the working environmentof a mobile network. Effectively the environment is modeled in twodimensions and rays are traced through the environment to determine thecoverage achieved. If the space around the environment is modeled as agrid, the process for drawing a line can be used to traverse a raythrough the grid. The line segment becomes a ray if no endpoint isassumed. As the line is traced through the environment, instead ofdrawing each pixel in the display bitmap, the grid cell can be checked.

Simulation of Sonar Propagation

In sonar propagation simulation, often laminar sheets of water or othermedia are traced with rays to determine the propagation of rays throughthe medium. If the media can be modeled on a discrete grid, the linedrawing algorithm can be used to propagate the ray through the medium.

Height Field Rendering

Ray tracing has long been used to render height fields or elevationmaps. The run-based ray line drawing technique can be adapted totraverse a ray over the height field to seek an intersection with anelevation.

Hough Transform

In line and shape detection problems in which the Hough Transform isused as a voting strategy, the voting process requires a line of votesto be cast into an array. The run-based line digitization techniquesdescribed can be used to cast the votes.

Additional Line Drawing Techniques that Use the Line StructureInformation

The line drawing techniques that are the subject matter ofPCT/US02/24711, the parent of this patent application, require thatinformation about the structure of the line being drawn be computed. Thestructural information can also be used to enhance other aspects ofdrawing a line. Three of these are discussed in the following:techniques for drawing lines or computing rays in parallel;techniques for adding anti-aliasing pixels to lines drawn using runs ofpixels; andtechniques for improving memory setting processes based on runs.Parallel Techniques for Determining the Cells of a Raster that areIntercepted by a Line: FIG. 16

An advantage of using runs of any order n to determine the cells of araster that are intercepted by a line is that once the first truncatedrun of order n is drawn and the structural information about the linehas been computed, it can be used to determine the positions in theraster of all of the other runs of order n of cells that are interceptedby the line. Consequently, if the hardware or software being used todetermine the cells permits it, the cells belonging to the other runs oforder n can be computed in parallel. The parallel computations can beperformed by calculating each run of order n in the cells intersected bythe line in parallel; or by calculating the positions of groups of runsof order n in the cells intersected by the line in parallel and usingthe sequential version of the algorithm to draw each group. Anotherlevel of parallelism could be attained within each run by settingmultiple runs of a lower order or pixels per iteration of the sequentialor parallel version of the algorithm. The run of order n could bedivided into runs of an order less than n, each of which could beprocessed in parallel.

FIG. 16 gives an overview of how runs may be computed in parallel. Thereare two parts: initialization 1615, which calculates and sets the firsttruncated run of cells and then computes the starting cells of each ofthe full runs of order i being set, and the part 1615 that sets the fullruns in parallel. Beginning with start 1603, the first part of theinitialization is done as shown in FIG. 6. As set forth in theexplanation of FIG. 6, truncated run processor 102 receives as inputsthe endpoints of the line; from the endpoints, it calculates the line'sslope and determines the starting pixel and the starting point and errorterm β for the first run truncated run of order 1. Then it proceeds asfollows for that order and all others 1 . . . n: first it calculates thelength of the truncated run of the order, then it sets the truncated runin the bitmap, and then it calculates the error term for the next order.At the end of the process, the location of the starting pixel of thefirst full run of order n is known, along with the error term. From thisinformation, the starting cells and structure parameters of each of thefull runs of order n may be computed, as shown at 1610. This completesthe initialization, and the full runs are set in parallel, with the lastrun terminating when the line's end point is reached.

The apparatus of FIG. 2 may be modified to perform the processing ofFIG. 16 by replacing bitmap processor [106] with a number of bitmapprocessors 106 (a . . . j), which can operate on bitmap memory 108 inparallel. Run processor 104 would then give each bitmap processor 106(i)the starting cell and structure parameters for the run that bitmapprocessor 106(i) is to set. The degree of order n may be determined bythe number of bitmap processors 106; for example, if there are fourbitmap processors, the degree of order n may be that at which the linecontains full runs and a final truncated or untruncated run such thatthe total number of the full and final runs is divisible by 4. Thedegree of order n may be determined by the designer or based on externalparameters such as the length and slope of the line. The bitmapprocessors may also be able to set an array of bits in a singleoperation; in that case, the order whose runs best fit the array thatthe bitmap processor can set might be the one for which the operation isparallelized.

Computing the Starting Points of the Runs for Order 1

To calculate the position of the j-th full length run in the line, letus reconsider the linear equation given byβ_(j) ^([1])−β₀ ^([1])=α^([1])(x _(j) ^([1]) −x ₀ ^([1]))−(y _(j) ^([1])−y ₀ ^([1])).where β₀ ^([1]) is the intercept value of the first full length rungiven at position (x₀ ^([1]), y₀ ^([1])). We know that j=y_(j) ^([1])−y₀^([1]) as the height of a run is always one.

Therefore the position of the j-th full length run in the line is(x _(j) ^([1]) , y _(j) ^([1]))=([β_(j) ^([1])−β₀ ^([1]) +j+x ₀^([1])]/α^([1]) , y ₀ ^([1]) +j)

To increase the efficiency of the algorithm, we can normalize each ofthe values at order 1 by α^([1]). This eliminates a number of divisionand multiplication operations.

Computing the Starting Points of the Runs for Order 2

To calculate the position of the j-th run in the line, let us reconsiderthe linear equation given byβ_(j) ^([1])−β₀ ^([1])=α^([1]() x _(j) ^([1]) −x _(0j) ^([1]))−(y _(j)^([1]) −y ₀ ^([1])).For order 2 we haveβ_(j) ^([2])−β_(o) ^([2])=α^([1]() x _(j) ^([2]) −x ₀ ^([2]))−α^([0]() y_(j) ^([2]) −y ₀ ^([2])).Define p_(j) ^([h]), where h≦i, to be the number of runs of order h inthe partial line up to the beginning of the j-th run of order i from thefirst full-length run of order 2. Thereforep _(j) ^([0])=x_(j) ^([2]) −x ₀ ^([2])p _(j) ^([1])=y_(j) ^([2]) −y ₀ ^([2])p_(j) ^([2])=jThe above linear equation for order 2 becomesβ_(j) ^([2])=α ^([1]) p _(j) ^([0])−α^([0]) p _(j) ^([1]).For runs of order 2, shape 0, we haveα^([2])=μ^([1])=(r^([1])+1)α^([1])−α^([0]) and by substituting forα^([0]) into the previous equation we haveβ_(j) ^([2])=α^([2]) p _(j) ^([1])−α^([1])((r ^([1])+1)p _(j) ^([1]) −p_(j) ^([0]))For this shape, each run of order 2 must contain a long run of order 1and r^([1])+1 short runs therefore α^([1])=(r^([1])+1)p_(j) ^([1])−p_(j)^([0]), andβ_(j) ^([2])=α^([2]) p _(j) ^([1])−α^([1]) p _(j) ^([2]).Similarly for shape 1 we getβ_(j) ^([2])=α_([1]) p _(j) ^([2])−α^([2]) p _(j) ^([1]).Hence for shape 0$y_{j}^{\lbrack 2\rbrack} = {y_{0}^{\lbrack 2\rbrack} + \frac{{j\quad\alpha^{\lbrack 1\rbrack}} - \beta_{j}^{\lbrack 2\rbrack} + \beta_{0}^{\lbrack 2\rbrack}}{\alpha^{\lbrack 2\rbrack}}}$$x_{j}^{\lbrack 2\rbrack} = {x_{0}^{\lbrack 2\rbrack} + \frac{y_{j}^{\lbrack 2\rbrack} + \beta_{j}^{\lbrack 2\rbrack} - \beta_{0}^{\lbrack 2\rbrack}}{\alpha^{\lbrack 1\rbrack}}}$and for shape 1$y_{j}^{\lbrack 2\rbrack} = {y_{0}^{\lbrack 2\rbrack} + \frac{{j\quad\alpha^{\lbrack 1\rbrack}} + \beta_{j}^{\lbrack 2\rbrack} - \beta_{0}^{\lbrack 2\rbrack}}{\alpha^{\lbrack 2\rbrack}}}$$x_{j}^{\lbrack 2\rbrack} = {x_{0}^{\lbrack 2\rbrack} + \frac{y_{j}^{\lbrack 2\rbrack} + \beta_{j}^{\lbrack 2\rbrack} - \beta_{0}^{\lbrack 2\rbrack}}{\alpha^{\lbrack 1\rbrack}}}$To increase the efficiency of the algorithm, we can normalize each ofthe values at order 2 by α^([2]), as in order 1. This eliminates anumber of division and multiplication operations. In this case we canalso define {circumflex over (α)}^([2])=min({circumflex over(v)}^([2]),{circumflex over (μ)}^([2])) to be the proportion of runs oforder i to order 1 where τ^([2]=)1/{circumflex over (α)}^([2]) is theaverage run length of order 2.Computing the Starting Points of the Runs for Order i

To describe a general algorithm for any order i we can expand thedefinition of the algorithm for order 2. The main differences betweenthe algorithm of order 2 and i, is that runs of any order i>2 occur inone of the four possible shapes not two. For orders beyond 2 howeverthere is no restriction that the first run in the line y=αx is long,therefore all four shapes are possible. Apart from this key difference,the calculation of the structural parameters, the length of the initialtruncated run length, the starting point of each run of order i and thelength of each run of order i can be derived from order 2.

To expand on our idea of shape at order 2, we define the type instead ofthe shape of the runs of order i>1 to be type t^([i])=0 ifμ^([i−1])<v^([i−1]), type t^([i])=1 if v^([i−1])<μ^([i−1]) and definethe type of order 1 to be t^([i])=0, we can calculate the shape by therelationships^([i])=t^([i−1])⊕t^([i])where the operator ⊕ denotes the binary concatenation of the two types.For example if t^([2])=1 and t^([3])=0, shape s^([3]) has the binaryvalue “10”, or decimal 2, and therefore, s^([3])=2: sl⁺. We can usethese definitions of type and shape for order 2 by assuming that thetype of all runs of order 1 will always be type 0.

The method of calculating the position of the j-th run of order i>2 alsofollows from our discussion of order 2. To calculate the position of thej-th run in the line, let us reconsider the linear equation given byβ₀ ^([1])−β₀ ^([1])α^([1])(x _(j) ^([1]) −x _(0j) ^([1]))−(y _(j) ^([1])−y _(0j) ^([1]))orβ_(j) ^([1])−β₀ ^([1])=α^([1]) p _(j) ^([0])−α^([0]) p _(j) ^([1]).

By induction, it is possible to show that for order i, type 0β_(j) ^([i])−β₀ ^([i])=α^([i]) p _(j) ^([i−1])−α^([i−1]) p _(j) ^([i])

-   and for order i, type 1    β_(j) ^([i])−β₀ ^([i])=α^([i−1]) p _(j) ^([i])−α^([i]) p _(j)    ^([i−1])

Therefore, because p_(j) ^([i])=j, we can substitute back to orders 1and 0 where p_(j) ^([1])=y_(j) ^([i]) −y ₀ ^([i]) and p_(j) ^([0])=x_(j)^([i]) −x ₀ ^([i]) in a similar manner as we performed for order 2.

To increase the efficiency of the algorithm, we can normalize each ofthe values at order i by α^([i]), similarly to what we did for order 2.This eliminates a number of division and multiplication operations. Inthis case we can also define {circumflex over (α)}^([i])=min({circumflexover (v)}^([i]), μ^([i])) to be the proportion of runs of order i toorder 1 where τ^([i])=1/α^([i]) is the average run length of order 2.

Implementing the Calculation of the Starting Positions and Compositionsof the Full Runs of Order 1

In the following, μ and v are the structural parameters described in theparent of the present application. For each scanline, j, in the imagegiven the line segment (x ₀, y₀) to (x ₁, y₁) such that └y₀┘≦j≦┌y₁┐where x₀, y₀, x₁, y₁ ∈ R and x₀<x₁:

-   1. If └y₀┘≦j≦┌y₁┐, the scan line intersects the line segment.-   2. The starting point for the line is:    x₀ ^([0])=└x₀┘    y₀ ^([0])=└y₀┘-   3. To handle the first pixel, the distance of the intercept to the    position of the first pixel is calculated:    β′=(y ₀ −└y ₀┘)−(x ₀ −└x ₀┘)α^([1])    If β′≧0 then β₀ ^([0])=β′. If β′21 0 then β₀ ^([1])=β′−α^([1]) and    the position of the first run is reset to (x ₀ ^([0]), y₀    ^([0]))=(└x₀+1, └y₀┘).-   4. Calculate the average run length of order 1.    $\tau^{\lbrack 1\rbrack} = \frac{x_{1} - x_{0}}{y_{1} - y_{0}}$-   5. The length of the initial truncated run is    r_(t)^([1]) = ⌈τ^([1]) − β̂₀^([0])⌉-   6. The intercept value at the start of the first full-length run of    order 1 is:    {circumflex over (β)}₀ ^([1])=r_(l) ^([0])−(τ^([1])−{circumflex over    (β)}β₀ ^([0]))-   7. If the scanline intersects the start point, i.e. if y₀ ^([0])=j,    set a run length r₀ ^([1]) at pixel (x ₀ ^([0]), y₀ ^([0])) and    cease.-   8. We now have └y₀┘<j≦┌y₁ ┐. Set the reference point for the line to    start of the first full-length run of order 1 in the line    (x ₀ ^([1]) , y ₀ ^([1]))=(x ₀ ^([0]) +r ₀ ^([1]) , y ₀ ^([0])+1)-   9. The length of a short run is    r_(s) ^([1])=└τ^([1])┘-   10. The portion of long runs of order 1 to the total number of run    of order 1 in the line is    {circumflex over (v)}^([1])=τ^([1]) −r _(s) ^([1])-   11. The vertical distance from the reference point to the intercept    of the line and the j-th run is    δ{circumflex over (β)}_(j) ^([1])={circumflex over (β)}₀ ^([1])    +j(1−{circumflex over (v)} ^([1]))-   12. The normalized intercept value of the line and the j-th run is    {circumflex over (β)}_(j) ^([1])=Δ{circumflex over (β)}_(j)    ^([1])−└Δ{circumflex over (β)}_(j) ^([1])┘-   13. The position of the run is    (x _(j) ^([1]) , y ^([1]))=({circumflex over (β)}_(j) ^([1]) +jτ    ^([1]) , j)-   14. If the scanline intersects the end point of the line segment,    i.e. if y₀ ^([1])=y₁, calculate the last truncated run length    r _(j) ^([1]) =y ₀ ^([1]) −┌y ₁┐    Set a run of length r_(j) ^([1]) at pixel (x _(j) ^([1]) , y _(j)    ^([1])) and cease.-   15. We now have └y₀┘<j<┌y₁┐ and only the intermediate runs remain.    The length of the j-th run is    r _(j) ^([1]) =r _(s) ^([1]) if {circumflex over    (β)}^([1])≧{circumflex over (v)}^([1])    r _(j) ^([1]) =r _(s) ^([1])+1 if {circumflex over (γ)} _(j)    ^([1])<{circumflex over (v)}^([1])    Set a run of length r_(j) ^([1]) at pixel (x ₁ ^([1]) , y _(j)    ^([1]))    Implementing the Calculation of the Starting Positions and    Compositions of the Full Runs of Order 2

We assume that the steps 1-5 and 7-10 have been executed and we havevalues for the parameters x₀ ^([0]) , y ₀ ^([0]), x₀ ^([1]) , y ₀^([1]), r_(t) ^([1]), r_(s) ^([1]), {circumflex over (β)}₀ ^([1]),{circumflex over (v)}^([1]), τ^([1]). If j=0 and the first possiblytruncated run in the line is to be set, we assume that the firstpossibly truncated run of order 1 has been set.

-   1. Calculate the shape of the run of order 1.    If 1−{circumflex over (v)}^([1])<{circumflex over (v)}^([1])    s^([2])=0    else    s^([2])=1-   2. Calculate the average run length of order 1.    If s^([2])=0    τ^([2])=1/(1−{circumflex over (v)} ^([1]))    else    τ^([2])=1/{circumflex over (v)} ^([1])-   3. The length of the initial truncated run is    If s^([2])=0    r _(t) ^([2])=|τ^([2])−{circumflex over (β)}₀ ^([1])|    else    r _(t) ^([2])=└{circumflex over (β)}₀ ^([1])┘-   4. The intercept value at the start of the first full-length run of    order 1 is:    If s^([2])=0    {circumflex over (β)}₀ ^([2]) =r _(t) ^([2])−(τ^([2])−{circumflex    over (β)}₀ ^([1]))    else    {circumflex over (β)}₀ ^([2])={circumflex over (β)}₀ ^([1]) −r _(t)    ^([2])-   5. If j=0, set a run of order 1 and length r₀ ^([1]) at pixel (x ₀    ^([0]) , y ₀ ^([0])). Then set a run of order 2, length r₀ ^([2])    and shape s^([2]) at pixel (x ₀ ^([1]) , y ₀ ^([1])) and cease.-   6. Set the reference point for the line to start of the first    full-length run of order 2 in the line:    If s^([2])=0    (x ₀ ^([2]) , y ₀ ^([2]))=(x ₀ ^([1])+r₀ ^([2])(r _(s) ^([1])+1)−1,    y ₀ ^([1])+r₀ ^([2]))    else    (x ₀ ^([2]) , y ₀ ^([2]))=(x ₀ ^([1])−r₀ ^([2])r_(s) ^([1])−1, y ₀    ^([1])−r₀ ^([2]))-   7. The length of a short run is    r_(s) ^([2])=└τ^([2])┘-   8. The portion of long runs of order 1 to the total number of run of    order 1 in the line is    If s^([2])=0    {circumflex over (v)}^([2])=τ^([2]) −r _(s) ^([2])    else    {circumflex over (v)}^([2])=1−(τ^([2]) −r _(s) ^([2]))-   9. The run-based distance from the reference point to the intercept    of the line and the j-th run is:    δ{circumflex over (β)}_(j) ^([2])={circumflex over (β)}₀ ^([2])    −j(1−{circumflex over (v)} ^([2]))-   10. The normalized intercept value of the line and the j-th run is    {circumflex over (β)}_(j) ^([2])=Δ{circumflex over (β)}_(j)    ^([2])−└Δ{circumflex over (β)}_(j) ^([2])┘-   11. The position of the run is    If s^([2])=0    y _(j) ^([2]) =y ₀ ^([2]) +jτ ^([2])+{circumflex over (β)}_(j)    ^([2])−{circumflex over (β)}₀ ^([2])    x _(j) ^([2]) =x ₀ ^([2]) +y _(j) ^([2])τ^([1])+({circumflex over    (β)}_(j) ^([2])−{circumflex over (β)}₀ ^([2]))/τ^([2])    else    y _(j) ^([2]) =y ₀ ^([2]) +jτ ^([2]){circumflex over (β)}_(j)    ^([2])+{circumflex over (β)}₀ ^([2])    x _(j) ^([2]) =x ₀ ^([2]) +y _(j) ^([2])τ^([1])−({circumflex over    (β)}_(j) ^([2])−{circumflex over (β)}₀ ^([2]))/τ^([2])-   12. The length of the j-th run is    r _(j) ^([2]) =r _(s) ² if {circumflex over (β)} _(j)    ^([2])≧{circumflex over (v)}^([2])    r _(j) ^([2]) =r _(s) ^([2])+1 if {circumflex over (β)}_(j) ^([2])    <{circumflex over (v)} ^([2])-   13. Set a run of length r_(j) ^([2]) shape s^([2])=t^([2]) at pixel    (x _(j) ^([2]) , y _(j) ^([2]))    Implementing the Calculation of the Starting Positions and    Compositions of the Full Runs of Order i

We assume we have values for the parameters x₀ ^([k]) , y ₀ ^([k]),r_(t) ^([k]), r_(s) ^([k]), {circumflex over (β)}₀ ^([k]), {circumflexover (v)}^([k]), τ^([k]), where k=0 . . . i−1. If j=0 and the firstpossibly truncated run of order i in the line is to be set, we assumethat the first possibly truncated run of order i−1 has been set.1. Calculate the shape of the run of order 1.If 1−{circumflex over (v)}^([i−1])<{circumflex over (v)}^([i−1])t^([i])=0elset^([i])=12. Calculate the average run length of order 1.If t^([i])=0τ^([i])=1/(1−{circumflex over (v)} ^([i−1]))elseτ^([i])=1/{circumflex over (v)}^([i−1])3. The length of the initial truncated run isIf t^([i])=0r _(t) ^([i])=|τ^([i])−{circumflex over (β)}₀ ^([i−[1])|

-   else    r_(t) ^([i])=└β₀ ^([i−1])┘-   4. The intercept value at the start of the first full-length run of    order 1 is:-   If t^([i])=0    {circumflex over (β)}₀ ^([i]) =r _(t) ^([i])−(τ^([i])−{circumflex    over (β)}₀ ^([i−1]))-   else    {circumflex over (β)}₀ ^([i])={circumflex over (β)}₀ ^([i−1]) −r    _(t) ^([i])-   5. If j=0, set a run of order i and length r₀ ^([i]) at pixel (x ₀    ^([i−1]) , y ₀ ^([i−1])).-   6. Set the reference point for the line to start of the first    full-length run of order 2 in the line:    δr^([k])=r₀ ^([i])-   for k=−1, . . . , 1    if t^([i])=0    δr ^([k])=δr^([k+1])(r _(s) ^([k])+1)−1    else    δr ^([k])=δr^([k]) r _(s) ^([k−1])−1    y₀ ^([1]) =y ₀ ^([1]) +δr ^([1])    x₀ ^([i]) =x ₀ ^([1]) +δr ^([0])-   7. The length of a short run is    r_(s) ^([i])=└τ^([i])┘-   8. The portion of long runs of order 1 to the total number of run of    order 1 in the line is    If t^([i])=0    {circumflex over (v)}^([i])=τ^([i]) −r _(s) ^([i])    else    {circumflex over (v)}^([i]) =1−(τ ^([i]) −r _(s) ^([i]))-   9. The run-based distance from the reference point to the intercept    of the line and the j-th run is:    δ{circumflex over (β)}_(j) ^([i])={circumflex over (β)}₀ ^([i])    +j(1−{circumflex over (v)} ^([i]))-   10. The normalized intercept value of the line and the j-th run is    {circumflex over (β)}_(j) ^([i])=Δ{circumflex over (β)}_(j)    ^([i])−└Δ{circumflex over (β)}_(j) ^([i])┘-   11. The position of the run is    p_(j) ^([1])=j    for k=i, . . . , 1    if t^([i])=0    p_(j) ^([k−1])=τ^([2]) p _(j) ^([k]) +(β _(j) ^([2])−β₀ ^([2]))    else    p_(j) ^([k−1])=τ^([2]) p _(j) ^([k])+(β_(j) ^([2])−β₀ ^([2]))    else    p_(j) ^([k−1])=τ^([2]) p _(j) ^([k])−(β_(j) ^([2])−β₀ ^([2]))    y_(j) ^([i]) =y ₀ ^([i]) −p _(j) ^([1])    x_(j) ^([i]) =x ₀ ^([i]) −p _(j) ^([0])-   12. The length of the j-th run is    r_(j) ^([i]) =r _(s) ^([i]) if {circumflex over (β)}_(j) ^([i])    ≧{circumflex over (v)} ^([i])    r_(j) ^([i]) =r _(s) ^([i])+1 if {circumflex over (β)}_(j) ^([i])    <{circumflex over (v)} ^([i])-   13. Set a run of length r_(j) ^([i]), shape s^([i])=t^([1−])⊕t^([i])    at pixel (x _(j) ^([i]) , y _(j) ^([i])).    Optimizations

For the parallel algorithm, we can rearrange the calculations to movemany of the expensive operations into the initialization phase of thealgorithm. For example Algorithm 1 illustrates a possible implementationof an order 1 algorithm.

Algorithm 1 Calculating the Length and Position of the j-th Run of Order1 in the Digital Line.

Constants pertaining to the line:$\tau^{\lbrack 1\rbrack} = \frac{1}{\alpha^{\lbrack 1\rbrack}}$c₀^([1]) = τ^([1])α^([0]) c₁^([1]) = τ^([1])β₀^([1])c₂^([1]) = τ^([1])v^([1])

Constants pertaining to the j-th run of order 1.$c_{3}^{\lbrack 1\rbrack} = {\frac{j}{\alpha^{\lbrack 1\rbrack}} = {j\quad c_{0}^{\lbrack 1\rbrack}}}$$c_{4}^{\lbrack 1\rbrack} = {\frac{\beta_{0}^{\lbrack 1\rbrack} + {j\quad\mu^{\lbrack 1\rbrack}}}{\alpha^{\lbrack 1\rbrack}} = {c_{1}^{\lbrack 1\rbrack} + {c_{3}^{\lbrack 1\rbrack}\mu^{\lbrack 1\rbrack}}}}$$c_{5}^{\lbrack 1\rbrack} = {\frac{\beta_{j}^{\lbrack 1\rbrack}}{\alpha^{\lbrack 1\rbrack}} = {c_{4}^{\lbrack 1\rbrack} - \left\lfloor c_{4}^{\lbrack 1\rbrack} \right\rfloor}}$

The position of the j-th run of order 1:x _(j) ^([1]) =x ₀ ^([1]) +c ₅ ^([1]) +c ₃ ^([1]) −c ₁ ^([1])y _(j) ^([1]) =y ₀ ^([1]) +j

The length of the j-th run of order 1:if c₅ ^([1])<C₂ ^([1])

-   -   -   The length of the j-th run of order 1 is long, r_(j)            ^([1])=|τ^([1])|.            else        -   The length of the j-th run of order 1 is short, r_(j)            ^([1])=└τ^([1])┘.

From Algorithm 1, the constants c₁ ^([1]) and c₂ ^([1]) can becalculated once in the initialization of the algorithm. The constants c₃^([1]) and c₄ ^([1]) can be calculated iteratively as each new paralleltask is created from preexisting constants. Therefore each run can beset in parallel for a cost of approximately 6 additions, amultiplication and a comparison compared to a cost of 3 additions andtwo comparisons for the central loop of the sequential algorithm. Hence,the overhead imposed by the parallelization of the algorithm is quicklyrecouped.

Similarly, for an order 2 implementation as described in Algorithm 1,migrating to a parallel implementation would add two divisions and amultiplication to the cost of initializing the algorithm and 10additions, a multiplication and division and a comparison to the cost ofcalculating the length and position of each run of order 2.

Algorithm 2 Calculating the Length and Position of the j-Th Run of Order2 in the Digital Line.

Constants pertaining to the line:$c_{0}^{\lbrack 2\rbrack} = \frac{\alpha^{\lbrack 1\rbrack}}{\alpha^{\lbrack 2\rbrack}}$$c_{1}^{\lbrack 2\rbrack} = \frac{\beta_{0}^{\lbrack 2\rbrack}}{\alpha^{\lbrack 2\rbrack}}$$c_{2}^{\lbrack 2\rbrack} = \frac{v^{\lbrack 2\rbrack}}{\alpha^{\lbrack 2\rbrack}}$

Constants pertaining to the j-th run of order 1:$c_{3}^{\lbrack 2\rbrack} = {{j\quad\frac{\alpha^{\lbrack 1\rbrack}}{\alpha^{\lbrack 2\rbrack}}} = {j*c_{0}^{\lbrack 2\rbrack}}}$$c_{4}^{\lbrack 2\rbrack} = {\frac{\beta_{0}^{\lbrack 2\rbrack} + {j\quad\mu^{\lbrack 2\rbrack}}}{\alpha^{\lbrack 2\rbrack}} = {c_{1}^{\lbrack 2\rbrack} + {c_{3}^{\lbrack 2\rbrack}\mu^{\lbrack 2\rbrack}}}}$$c_{5}^{\lbrack 2\rbrack} = {\frac{\beta_{j}^{\lbrack 2\rbrack}}{\alpha^{\lbrack 2\rbrack}} = {c_{4}^{\lbrack 2\rbrack} - \left\lfloor c_{4}^{\lbrack 2\rbrack} \right\rfloor}}$if μ^([1])<v^([1])

-   -   -   The length of the j-th run of order 2, shape 0: l⁺s:        -   if c₅ ^([2])<c₂ ^([2])            -   The length of the j-th run of order 2 is long, r_(j)                ^([2])=|τ^([2])|.        -   else            -   The length of the j-th run of order 2 is short, r_(j)                ^([2])=└τ^([2])|.                The height of the j full runs of order 2, shape 0: l⁺s:                y=c₅ ^([2])+c₃ ^([2])−c₁ ^([2])                else        -   The length of the j-th run of order 2, shape 1: ls⁺:        -   if c₅ ^([2])<c₂ ^([2])            -   The length of the j-th run of order 2 is short, r_(j)                ^([2])=└τ^([2])┘.        -   else            -   The length of the j-th run of order 2 is long, r_(j)                ^([2])=|τ^([2])|.        -   The height of the j full runs of order 2, shape 1: ls⁺:            y=c₅ ^([2])−c₃ ^([2])+c₁ ^([2])            The position of the j-th run of order 2 (x _(j) ^([2]) , y            _(j) ^([2]))            $x_{j}^{\lbrack 2\rbrack} = {x_{0}^{\lbrack 2\rbrack} + \frac{y + \beta_{j}^{\lbrack 2\rbrack} - \beta_{0}^{\lbrack 2\rbrack}}{\alpha^{\lbrack 1\rbrack}}}$            y_(j)^([2]) = y₀^([2]) + y            Algorithm 3 Calculating the Length and Position of The j-th            Run of Order i in the Digital line.

Constants pertaining to the line:$c_{0}^{\lbrack i\rbrack} = \frac{\alpha^{\lbrack{i - 1}\rbrack}}{\alpha^{\lbrack i\rbrack}}$$c_{1}^{\lbrack i\rbrack} = \frac{\beta_{0}^{\lbrack i\rbrack}}{\alpha^{\lbrack i\rbrack}}$$c_{2}^{\lbrack i\rbrack} = \frac{v\lbrack i\rbrack}{\alpha^{\lbrack i\rbrack}}$

Constants pertaining to the j-th run of order i:$c_{3}^{\lbrack i\rbrack} = {\frac{j}{\alpha^{\lbrack i\rbrack}} = {j\quad c_{0}^{\lbrack i\rbrack}}}$$c_{4}^{\lbrack i\rbrack} = {\frac{\beta_{0}^{\lbrack i\rbrack} + {j\quad\mu^{\lbrack i\rbrack}}}{\alpha^{\lbrack i\rbrack}} = {c_{1}^{\lbrack i\rbrack} + {c_{3}^{\lbrack i\rbrack}\mu^{\lbrack i\rbrack}}}}$$c_{5}^{\lbrack i\rbrack} = {\frac{\beta_{j}^{\lbrack i\rbrack}}{\alpha^{\lbrack i\rbrack}} = {c_{4}^{\lbrack i\rbrack} - \left\lfloor c_{4}^{\lbrack i\rbrack} \right\rfloor}}$

The shape of the run of order i:if μ^([i−1])<v^([i−1])t^([i])=1elset^([i])=1

The shape of the run of order i:s^([i])=t^([i−1])⊕t^([i])

The length of the j-th run of order i:if t^([i])=0if c₅ ^([i])<c₂ ^([i])

-   -   -   -   The length of the j-th run of order i is long, r_(j)                ^([i])=|τ^([i])|.

        -   else            -   The length of the j-th run of order i is short, r_(j)                ^([i])=└τ^([i])┘.                else

        -   if c₅ ^([i])<c₂ ^([i])            -   The length of the j-th run of order i is short, r_(j)                ^([i])=└τ^([i])┘.

        -   else            -   The length of the j-th run of order i is long, r_(j)                ^([i])=|τ^([i])|.

The start point of the j-th run of order i:p_(j) ^([i])=jfor k=i−1 to 1if t^([k])=0p _(j) ^([k−1])=α^([k−1]) p _(j) ^([k])/α^([k])+(β_(j)^([k])−β^([k]))/α^([k])=c₀ ^([k]) p _(j) ^([k]) +c ₅ ^([k]) −c ₁ ^([k])elsep _(j) ^([k−1])=α^([k−1]) p _(j) ^([k])/α^([k])−(β_(j) ^([k])−β₀^([k]))/α^([k]) =c ₀ ^([k]) p _(j) ^([k]) −c ₅ ^([k]) +c ₁ ^([k])(x _(j) ^([i]) , y _(j) ^([i]))=(p _(j) ^([0]) +x ₀ ^([i]) , p _(j)^([1]) −y ₀ ^([i]))Efficiency of Parallelization

To compare the efficiency of the parallel and sequential versions of thealgorithms, we carry out an asymptotic analysis of the time complexityof the overall hierarchy algorithm for runs to all orders. At the outsetwe point out that we assume a rational slope α=a/d . We assume a and dare positive integers with a<d. For the case of an algebraic irrationalslope α as described in the parent of the present application, adiscussion of computing primitives would be required. However, givensuch a discussion, most of the rational slope analysis can be carriedover to the algebraic irrational case. This follows because all of theoperations involved in the algorithm are rational and it is well knownthat equations and inequalities in algebraic irrational numbersinvolving explicitly rational (in fact algebraic) operations can bedecided.

It has been shown in P. D. Stephenson and B. Litow, “Why step when youcan run: Iterative line digitization algorithms based on hierarchies ofruns”, IEEE Computer Graphics and Applications, 20(6):76-84,November/December 2000) that the number of orders in the hierarchy ofruns for a line with slope α is bounded above by O(log d) and that thenumber of runs in successive orders diminishes very much like partialquotients in the GCD algorithm for a and d. We will make use of thisimportant observation in our analysis.

Initializations across orders are sequential, i.e., values for order idepend on values for order i−1. It is evident that at each order,initialization costs O(1) rational arithmetic operations. We emphasizethat the division occurs only within the scope of the floor function.This suggests that the precision needed for division can be reducedconsiderably below the obvious requirement.

Within each order the computation of run lengths is actually aninitialization and we count it as such. The run positions can becomputed in parallel. Each position costs O(1) rational arithmeticoperations. However, here division appears to require full precision.Assembling our accounting, we have that initialization over all orderscosts O(log d) rational arithmetic operations. Idealizing slightly thepartial quotients behavior of GCD, we credit order i with d/2^(i) runs.If we allocate d/2^(i) processors to the i-th run, the time to computerun positions across all orders is O(1) rational arithmetic operations.This gives a work (time×processors) value of O(d)·(1+1/2+ . . .+1/2^(h)), where h is the maximum order, h≦log d. This yields a workvalue of O(d). If we allocate fewer processors, the parallel timeincreases, but it should be possible to choose configurations thatpreserve the O(d) work value.

In computer graphics and even computational geometry applications of thealgorithm, the value of d will probably remain below, say, 10⁴. Theboundary condition here is the realization of the digitized line. Ingraphics the constraint that d<10 ⁴ seems realistic, since we would beresolving beyond 10⁸ pixels were d to exceed 10⁴. For computationalgeometry applications, more thought about implicit representations ofdigitized objects is needed. Using hierarchical run-based methods hasshown a significant improvement in the cost to digitize a line. Animportant reason for the improvement has been the structural informationof the digital line provided by its run-based description.

To implement a system based on one of the algorithms presented, theorder of the algorithm must be decided either by the designer ordynamically by the system itself. While the length of the runs in thedigital line increases exponentially with the order of the run, the costof initializing an algorithm also rises, albeit linearly. Given that thecost of an iteration to decide the position and length of the next runin the line is so small, the number of runs in the line must besignificant enough to warrant using a higher order algorithm.

Problems of Ray Traversal in Three Dimensions: FIG. 17

There are many situations where a program needs to compute the traversalof a ray through three-dimensional space. One example of such a programis one that produces a two-dimensional display of a three-dimensionalobject; another is one that determines how a signal is propagatedthrough a volume. In both cases, there are two problems: representingthe traversal in the program and accelerating the traversal throughthose parts of the three-dimensional space that do not affect the ray. Agood general discussion of ray traversal in three dimensions can befound in Marc Levoy. Efficient ray tracing of volume data. ACMTransactions on Graphics, 9(3): 245-261, July 1990. FIG. 17, which showshow ray traversal may be used to make a two-dimensional representationof a three-dimensional object, is a slightly modified version of Levoy'sFIG. 2.

The two-dimensional representation of the three-dimensional object is animage 1705 that contains P×P pixels. The image is effectively a “window”into an image space 1701 that contains P×P×W voxels. A voxel is thethree-dimensional equivalent of a pixel: it is a data item thatcorresponds to a cell located at a point (x,y,z) in image space. Thedata in the voxel contains information about the cell; an example ofsuch information is the color and opacity of the cell. If the color andopacity correspond to data values, the data in the voxel may alsoinclude the data values. If image space 1701 is non-empty, it willinclude one or more object spaces 1703. Each object space has as manyvoxels as are required to contain all of the voxels for a particularobject in image space 1701. Each voxel 1709 in object space 1703contains data as required for the portion of the object represented bythe voxel. It should be noted that, as shown in FIG. 2, object space1703 may not have the same x, y, and z-axes as image space 1701, andthat consequently, a ray 1711 may pass through object space 1703 at anyangle relative to the axes of object space 1703.

To make image 1705, one or more rays 1711 corresponding to one or morepixels 1715 in the image are traced through image space 1701. If ray1711 passes through a voxel 1709 whose data affects the appearance ofthat pixel, the relevant information from the voxel becomes part ofsample 1713 for the pixel. When ray 1711 has completely traversed imagespace 1701, the information in sample 1713 is used to determine theappearance of pixel 1715.

As can be seen from FIG. 17, an image space 1701 of any size has anenormous number of voxels, and most of the voxels traversed by aparticular ray 1711 are either not part of an object space 1703, andtherefore provide no information to the traversing ray, or are behind anopaque voxel and can provide no information to the traversing ray beyondwhat was obtained from the opaque voxel. Efficient traversal of imagespace 1701 by a particular ray thus requires both efficient techniquesfor computing the voxels that are intersected by the ray and efficienttechniques for determining which of the voxels intersected by the raycontain information that is relevant to the pixel that corresponds tothe ray. As will be explained in detail in the following, the run-basedtechniques described in the parents of the present patent applicationmay be employed both to provide efficient computation of the voxels thatare intersected by a ray and to provide efficient determination of whichof the voxels intersected by the ray contain information relevant to theray's pixel. The same techniques can of course be employed in otherapplications where the information gathered By the ray traversal hasuses other than or in addition to display generation.

Using Run-Based Techniques to Compute Voxels Intersected by a Ray: FIGS.18-22

Geometry of Ray Projection: FIG. 18

The techniques that are explained in detail in the following are basedon the fact that a ray which is traversing a volume in which positionsof points are defined relative to X, Y, and Z axes may be represented astwo projections: one on the plane defined by the X and Y axes andanother on the plane defined by the X and Z axes. This is shown in thediagram of FIG. 18.

FIG. 18 shows a volume 1801 that is being traversed by a ray 1816.Points in volume 1801 are specified by the triplet (x,y,z), where x is avalue on X axis 1803, y is a value on Y axis 1805, and z is a value on Zaxis 1807. X axis 1803 is termed the major axis with regard totraversing ray 1816 because the direction of traversing ray 1816 iscloser to that of the X axis 1803 than it is to the other two axes.Expressed mathematically, X axis 1803 is the major axis because for anytwo points on ray 1816, (a,b,c) and (i,j,k), the difference δ between aand i will be greater than the difference δ between b and j or c and k.

Traversing ray 1816 passes through the top and bottom of volume 1801.The top is defined by the y,z plane 1815 made up of all of the points(x,y,z) where x=x(c). The bottom is defined by the y,z plane 1813 ofpoints ( x,y,z) where x=x(0). The point at which ray 1816 intersectsplane 1815 is point (x(c), y(d), z(e)) 1819; the point at which the rayintersects plane 1813 is ( x(0),y(a), z(b)) 1878. Two other planes ofinterest are x,y plane 1833 made up of all points for which z=z(0) andx,z plane 1835 made up of all points for which y=y(0). These planes areperpendicular to top plane 1815 and bottom plane 1813 and intersect at Xaxis 1803. As shown in FIG. 18, ray 1816 may be projected onto planes1833 and 1835. One way of doing the projection onto plane 1833 is toconstruct a line on top plane 1815 which passes through point 1819 andis perpendicular to the edge of plane 1833 which intersects plane 1815,construct a line on bottom plane 1813 which passes through point 1817and is perpendicular to the edge of plane 1933 which intersects plane1813, and then construct the projection 1821 of the ray on plane 1833 byconnecting the points at which the perpendicular lines intersect theedges. Projection 1827 on plane 1835 can be constructed in the samefashion.

Between them, projections 1821 and 1827 contain all of the informationneeded to determine the (x,y,z) coordinates of all of the points involume 1801 that are on ray 1816; consequently, projections 1821 and1827 together may be used to represent ray 1819. In the same fashion, athree-dimensional lattice of voxels in volume 1801 may be represented bytwo-dimensional lattices in planes 1833 and 1835 and the voxelsintersected by ray 1816 may be determined from the lattice cellsintersected in plane 1833 by projection 1821 and the lattice cellsintersected in plane 1835 by projection 1827.

Using Run-Based Techniques to Compute Voxels that are Intersected by aRay: FIGS. 19-22

Determining cells in a two-dimensional lattice that are intersected by aline is of course the problem addressed by the techniques described inthe portion of the present patent application that are carried over fromthe present application's parents. As set forth there, given two pointswithin a two-dimensional lattice, it is possible to describe latticecells that are intersected by a line connecting the two points in termsof runs and runs of runs of cells. If the lattices are projection planessuch as planes 1833 and 1835 and the lines are projections of a ray onthe projection planes, then it is further possible to determine voxelsintersected by the ray from the lattice cells intersected by theprojection lines.

How voxels intersected by the ray are determined depends upon whetherthe intersected voxels may be related to each other as 26-connectedvoxels or must be related to each other as 6-connected voxels.Six-connected voxels must connected to each other face-to-face; since agiven voxel has six faces, and it may be connected face-to-face with atmost six other voxels; hence the term “six-connected”. 26-connectedvoxels may also be connected to each other by any of their 8 corners and12 edges, hence the term “26-connected”. For volume rendering, theminimal 26-connected path of a ray suffices for many applications. Incomparison, iso-surfacing volumetric data and ray tracing require the6-connected ray path. In ray tracing, including edge and cornerintersections can also be beneficial to avoid possibly expensiveray-object intersection tests.

Determining 26-Connected Voxels Intersected by a Ray: FIGS. 19 and 20

FIG. 19 shows how to determine 26-connected voxels intersected by a rayfrom runs of lattice cells on the XY plane 1833 and the XZ plane 1835that are intersected by the projections of the ray on the XY and XZplanes. Volume 1801 in FIG. 19 has an X-axis 1803, a Y-axis 1805, and aZ-axis 1807, with the X-axis being the major axis. Projections of theray have been made on XY plane 1833 and on XZ plane 1835, and thetechniques described above have been used to determine lattice cells1927 intercepted by the projection on the XY plane and lattice cells1921 intercepted by the projection on the XZ plane. Lattice cells 1921and 1923 are then used to determine voxels 1909 intersected by the ray.

Determining the runs in the XY plane and the XZ plane is doneiteratively or in parallel a run at a time with the corresponding runsof each traversal determining the length of the resultant run of voxels.Starting at the point where the ray intersects bottom plane 1813, theinitial runs of lattice cells intercepted by the projections onto the XYand XZ planes are calculated. The shorter of these two runs and theequivalent number of cells of the longer run map to the first run ofvoxels intercepted by the ray. Thus, in FIG. 19, run of lattice cells1923 and the corresponding portion of run 1929 map onto run of voxelcells 1911. The remaining portion of the longer run 1929 and acorresponding portion of run 1925 map onto run of voxel cells 1913. Theremaining portion of run 1925 and all of run 1931 are then mapped ontorun of voxel cells 1915. Voxels intersected by the three-dimensionalline can therefore be determined by taking the length of the shorter ofthe two corresponding runs of cells as the length of thethree-dimensional run of voxels. The longer two-dimensional run of cellsis then reduced by the length of the shorter run, which is replaced bythe next two-dimensional run in its traversal. The procedure is moreprecisely described in Algorithm 2001 in FIG. 20. The cells in thetwo-dimensional runs need only be 8-connected, i.e., by any of theirfour sides or four corners.

Determining the 6-Connected Voxels Intersected by a Ray: FIG. 21

To construct the 6-connected path of the ray from the ray's projectionson primary and secondary projection planes 1833 and 1835, all that isneeded is to include an extra cell at the beginning of eachtwo-dimensional run. This process is described in FIG. 21. The runs 1919and 1931 in plane 1833 and the runs 1923 and 1925 are calculated asdescribed for FIG. 19; however, in set of runs 2101, a cell 2107 isadded at the corner where run 1929 and 1931 meet, and the same is donein set of runs 2109.

The 6-connected set of voxels that the ray intersects is shown at 2117.Run 2119 corresponds to run 1923 and the corresponding portion of run1929; run 2121 corresponds to the remaining portion of run 1929, addedcell 2115, and the part of run 1925 corresponding to the remainingportion of run 1929. Run 2123 corresponds to the remaining portion ofrun 1925, added cell 2107, and run 1931. All of the voxels in set ofvoxels 2117 are connected to one or more other voxels in the set bytheir faces.

Edge and Corner Intersections in Sets of 6-Connected Voxels: FIGS. 22-25

With 6-connected sets of voxels that are intersected by a ray,ray-object intersection tests can be avoided if voxels in which the rayintersects only an edge or a corner of the voxel can be excluded fromthe 6-connected set of voxels. FIG. 22 shows the edge intersection case.First, the edge intersection case must be detected. For an edgeintersection to occur in the set of 6-connected voxels intercepted bythe ray, a corner intersection must occur in one of the two-dimensionalsets of runs. This is shown in FIG. 22, where two-dimensional set ofruns 2205 that intersects projection 2201 has a corner intersection 2209between runs 2207 and 2211, but two-dimensional set of runs 2213 thatintersects projection 2203 does not have. In terms of the projections2201 and 2203, the value of β at point 2209 is 0, while its value at2218 is greater than 0. The edge intersection appears at 2225. Whenthere is an edge intersection, the voxel corresponding to cell 2220 isnot included in set of voxels 2219.

To detect an edge intersection, one proceeds as follows: when making thenext run of cells intersected by the projection, the distance 3 betweenthe projection and the corner of the next run is considered. If thedistance is non-zero, no edge intersection occurs. If it is zero, anedge intersection occurs and no voxel corresponding to the cell that isadded to the beginning of the next two-dimensional run is added to theset of voxels. The algorithm is set forth more precisely at 2301 in FIG.23, where {circumflex over (β)} indicates a normalized value of β

For a corner intersection to occur in the set of voxels, there must bean edge intersection at the same place in the sets of cells thatintersect the projections in planes 1833 and 1835. Therefore, a cornerintersection in the set of voxels can only occur when the two projectedrun lengths are the same length. This is shown in FIG. 24. Projectedlines 2401 and 2403 produce sets of two-dimensional runs 2405 and 2413in which runs 2407 and 2415 are the same length and for which δ for thenext run is 0. This results in corner intersections 2409 and 2417. Setof voxels 2421 then has a run of voxels 2423 that corresponds totwo-dimensional runs 2407 and 2415 and a run of voxels 2427 thatcorresponds to two-dimensional runs 2411 and 2419. Runs of voxels 2427and 2423 are connected at corner 2425.

To detect runs of voxels that are connected at a corner, one proceeds asfollows: If no edge intersection occurs, to get from the end of theprevious run of voxels (x′, y, z) to the start of the next run ( x′,y+1, z+1), the ray must intersect the voxel (x′, y+1, z) or (x′, y, z+1)as shown in FIG. 22. Which cell can be determined from the relativevalue of the normalized intercept values {circumflex over (β)} asdescribed in algorithm 2501 of FIG. 25. In this algorithm, the positionof the run has already been updated to the start of the next run.Therefore, the voxel intersected by the ray must be handled inretrospect once a corner intersection has been identified.

Efficiency of Computing Voxels Intersected by a Line fromTwo-Dimensional Runs of Cells Intersected by Projections of the Line

The three-dimensional traversal algorithm requires two separatetwo-dimensional traversals of each of the projections of the lineperpendicular to the major axis. If corner and edge intersections areignored, the number of runs in the three-dimensional path of the ray isthe sum of the number of runs in the two two-dimensional traversals.Therefore the order of the three-dimensional ray traversal algorithm isO(n/τ), where τ is now the average length in the three-dimensional pathof the ray and n is the associated number of cells.

Using Descriptions of Rays as Sets of Runs of Voxels to Accelerate RayTraversal: FIGS. 26-31

Using Encoding Runs to Describe the Voxels in a Volume: FIG. 26

A problem in rendering the contents of a volume as pixels by traversingthe volume with rays to obtain the information needed to define theappearance of the pixels is that in general, not all contained in thevolume along the path traversed by a particular ray has anything to dowith the appearance of the pixel. Two examples of this fact are thefollowing:voxels that represent empty space can have no effect on the appearanceof the pixel.once the ray that defines the appearance of the pixel has encountered avoxel that represents an opaque surface, the values of the furthervoxels traversed by the ray will have no effect on the appearance of thepixel.

The second idea, called early ray termination, only applies to certainrendering functions such as shaded rendering and isosurfaces. If theentire path of the ray through the volume must be determined, as for anX ray style or maximum intensity projection (MIP) rendering, early raytermination is not valid.

A general solution to the problem of skipping uninteresting space is tofind some way of indicating that the content of a particular set ofvoxels will have no effect on a traversing ray. The technique generallyused is to subdivide the space being traversed into regions that willaffect a traversing ray and regions that will not. When the ray enters aregion that does not affect it, there is no need to examine individualvoxels in the region.

The technique used to subdivide the space into regions in a preferredembodiment is similar to the compression technique termed run-lengthencoding. Run-length encoding is used in situations where the data beingcompressed includes sequences of data in which the same value isrepeated over and over again. Such a sequence of data may then berepresented by giving the value and the number of repetitions. Forexample, the sequence of 10 eights, 8888888888 may be encoded as thepair {value, number of occurrences}, or {8,10}. The compressed data thusbecomes a list of pairs {value, number of occurrences}. When thesequence of data represents pixels in an image, there may be a list ofpairs corresponding to each row of pixels in the image. A row of pixelsis of course a run, and in the following, runs used for encoding will betermed encoding runs, to distinguish them from the runs of cellstraversed by a line or ray.

Similarly, when the sequence of data represents voxels in a volume, theencoding runs may be rows of voxels in the volume. FIG. 26 gives anexample of an encoding run 2603 in volume 1801. Encoding run 2603 is arow of voxels in 1801 for which the coordinates are (x(i),y(0 . . .n),z(j)). Similar encoding runs may of course be made where the voxelsall have the same y and z coordinates and a range of x coordinates orthe same x and y coordinates and a range of z coordinates. The bestencoding runs to use with a particular ray are those for which thedimension with the range of coordinates is the dimension of the majoraxis of the ray. In encoding run 2603, there are only two areas 2605 and2607 in which the contents of the voxels will affect a ray thattraverses the voxel. Encoding run 2603 can thus be encoded as a list ofsets of three values as follows: {NULL,y(0),y(a)}; {VAL1, (y(b),y(b)};{NULL,y(c),y(d)}; {VAL2, (y(e),y(f)}; {NULL,y(g),y(h)}, where a part ofthe run that has no effect on the ray is indicated by NULL followed bythe y coordinates of the beginning and end of the sequence of voxelsthat have no effect and a part of the run that can affect a ray isindicated in the same way, except that NULL is replaced by value(s) thatwill affect the ray in the sequence.

Using Runs of Voxels Traversed by a Ray Together with Encoding Runs toAccelerate Ray Traversal: FIGS. 27-30

If lists of encoding runs 2603 have been made for each of the principleaxes of a volume, then the problem of determining how a ray is affectedby the voxels of the volume can be solved by determining whether thereare any voxels traversed by the ray (termed in the following ray voxels)in the list of encoding runs for the ray's major axis that are the sameas voxels which the list of encoding runs indicate will affect the ray.These later voxels belong to what will be termed in the followingsignificant runs of the encoding run. If there are ray voxels that arethe same as voxels in a significant run, the data items corresponding tothe voxels in the significant run that match the ray voxels must beexamined to determine how they affect the ray. Using the techniques justexplained to determine the ray voxels for a ray make determining whetherthere are ray voxels that are the same as voxels in a significant run ina particular encoding run easier, since in order for there to be suchvoxels, a run of the ray voxels must intersect the significant run.Further, as already described, the projections of the ray may bedescribed using runs of orders greater than 1, and the runs of rayvoxels may be thus made of any useful length. In general, of course, theemptier the volume is, the more acceleration results from the use oflonger runs of ray voxels.

How the ray voxels intersect significant runs is shown in FIG. 27. FIG.27 shows a slice 2701 of a volume of voxels in which the voxels areassociated with data from the Visible Human Project male dataset. Thedata describes the head and neck of a human male. Slice 2701 is a set ofvoxels that have the full range of x and y coordinates permitted by thevolume and a single z coordinate. The interesting region of slice 2701is represented as a list of encoding runs. The run of voxels representedby one such run, 2703(1), is shown in FIG. 27. As shown there, run2703(i) has three runs 2707 that contain no information that will affecta ray and two significant runs 2705; the encoding run corresponding torun of voxels 2703(i) will represent that fact as described above. Asingle ray 2709 is shown traversing the volume. The ray has the X-axisas its major axis; consequently, the encoding runs of interest are thosefor the X-axis. The ray is further parallel to the z-axis, and is thuscompletely contained in the encoding runs of slice 2701 throughout itstraverse of the volume. The techniques described here can also be usedin the case where the ray intersects more than one slice. Ray 2709intersects a set of runs of ray voxels; one of the runs of ray voxels isshown at 2711; this run also intersects an encoding run 2703(1) whosevoxels contain information that will affect ray 2709.

The problem of finding an intersection between ray 2709 and asignificant run 2705 of an encoding run 2703 can be solved by detectingany intersection between the run from the ray and a significant run inthe list of runs associated with the encoding run. Given that the listof significant runs from the encoding run is ordered, the search can beperformed efficiently and a number of strategies are possible. Thesimplest is described in algorithm 2801 of FIG. 28. For each run of rayvoxels belonging to the ray, if the end of the run of ray voxels isbefore the start of the significant run being checked, there is no rayvoxel in the run which is identical with a voxel of the significant runand a failing result can be returned for that run of ray voxels. On theother hand, if the start of the run of ray voxels is before the end ofthe significant run being considered, an intersection has occurred andthe voxels of the significant run that also belong to the run of rayvoxels must be considered in detail. In algorithm 2801, a linear searchis used to locate the intersecting ray voxels and significant runvoxels. As shown in algorithm 2901 in FIG. 29, once it has beendetermined that there is an intersection, binary search techniques canbe used to determine the intersecting ray voxels and significant runvoxels

Using Statistical Information about a Significant Portion of an Encodingrun to Accelerate a Ray Traversal

Further acceleration of ray traversal may be achieved by hierarchicallygrouping encoding runs that lie in the same or adjacent beams of thevolume data. A beam is the maximum set of voxels in a slice with thesame abscissa or ordinate. Therefore, voxels are arranged into beams,beams are arranged into slices, and slices are arranged into the volume.If the encoding run has already been divided into significant runs,aggregate information may be maintained about each significant run. Forexample, for a maximum intensity projection (MIP) rendering of avolumetric data set, storing the maximum value found within eachsignificant run would be useful. For a MIP rendering, as the rayintersects a significant run, the maximum value of the run can becompared against the current maximum for the ray. If the ray is alreadyat or above the maximum intensity, there is no need to consider thevoxels in the significant run in detail. In many cases, including suchstatistical information about the run dramatically reduces the number oftimes the original data set has to be accessed. The use of aggregateinformation as just described establishes a hierarchy of informationconcerning the encoding run and permits each partition of an encodingrun to be checked in its entirety for relevance to the ray before theintersecting voxels of the ray and the significant run are checked indetail.

Another example is isosurface rendering in which the surface thatcorresponds to a given value is reconstructed from the volumetric data.To reconstruct an isosurface using ray-tracing, each ray is tracedthrough the volume and where it would cross the isosurface iscalculated. Therefore if the viewer is external to the volume and theexterior of the volume is assumed to be empty space, the voxel in whichone of its neighboring values is greater than the isosurface value issought. Within this voxel, the ray may cross the isosurface, and if so,the point of the intersection must be calculated. To quickly skip overregions of the volume in which no isosurface intersection can exist, theminimum and maximum values found within each significant run can bekept. These values can be used to determine if an isosurfaceintersection can exist in any composite voxel. If not, the entiresignificant can be skipped.

Efficiency of Run-Length Encoding of Voxels

A three-dimensional run length encoding of a volume provides an O(kn²)compression of the volumetric data where k is proportional to the meannumber of significant runs in the encoding run lists and n³ is the sizeof the volume. Such runs are termed in the following significant runs.The maximum number of significant runs in list cannot exceed n/2 bydefinition. To traverse the run length encoding of the volume using therun-based ray traversal scheme, each run in the path of the ray isintersected with the significant runs in the encoded run list. As it wasimplemented, the searching strategy is an O(k) operation. The procedurein its entirety is therefore O(n/τ+k). Hence at its worst the behaviorof the algorithm is O(n) when the absolute value of the slope of the twotwo-dimensional projections of the ray is 1. Typically the behavior willbe O(n/τ) as k will typically be much less than n/2.

Organizing an Encoding run into Sets of Partitions of Various Sizes:FIGS. 30 and 31

The three-dimensional ray traversal algorithm described is veryefficient and for most of the rendering functions employed, such asisosurfacing and maximum intensity projection, the cost to process eachvoxel traversed by the ray is not high. The aim of using spatialsubdivisions to accelerate the traversal process is to reduce the amountof the volume that has to be processed and to concentrate the processingperformed to that part of the volume that affects the final image.However, if traversing the ray through the spatial subdivisionintroduces too much overhead, its contribution is limited.

Since only encoding voxels that are intersected by a run of ray voxelsbelonging to a particular ray are of interest for the ray's traversal,it is useful to divide the encoding runs into sets of partitions ofdifferent sizes. Aggregate information for each partition is maintained,and when a ray traverses a volume, the size of partition that best fitsthe orientation of the ray is employed.

This is shown in FIG. 30. Shown at 3001 is again a slice of a volume ofvoxels in which the voxels are associated with data from the VisibleHuman Project male data set. The data describes the head and neck of ahuman male. Slice 3001 is a set of voxels that have the fill range of xand y coordinates permitted by the volume and a single z coordinate. Ray2709 intersects slice 3001. The voxels of ray 2709 occur in runs of twolengths: a short run of 3 voxels is shown at 2711 and a long run of fourvoxels is shown at 2612. The encoding runs 3003 have been divided intorun length partitions. The partitions 3005 at the lowest level of thehierarchy each have 4 voxels; each higher-level partition in thehierarchy would contain a greater number of voxels than the previouslevel. In FIG. 30, the increase in size is based on the powers of two. Apartition at a given level k of the hierarchy contains 2^(k) of thelowest-level partitions. Thus, in FIG. 30, partition 3005 contains fourvoxels, partition 3007 contains two partitions 3005 and therefore 8voxels, and so forth for 1≦k≦4. While the partitions of the hierarchycan be related to each other in other ways, use of the powers of 2 makesit easy to calculate the partition sizes.

Ray 2709 has runs of 3 and 4 voxels; consequently, the best size ofpartition to use for ray 2709's traversal of slice 3001 is partition3005. Once the proper partition size has been selected, the algorithm isas shown at 3101 in FIG. 31. There are two cases: if the portion of theray run being processed extends beyond the partition and if it does not.In the first case, the partition's aggregate data is consulted todetermine whether the partition is relevant to the ray; if it is, thevoxels of interest in the partition are traversed. Thereupon, the lengthof the ray run is decreased by the number of intersecting voxels, andthe partition indicator is set to the next partition. In the secondcase, the partition's aggregate data is again consulted to determinewhether the partition is relevant; if it is, the voxels of interest inthe partition are traversed. Thereupon, the amount of the partition tobe processed is reduced by the length of the ray run.

Given the run-based ray traversal algorithm, the problem is thereforehow best to traverse the ray through the run length partitioning andwhat is the best partition size to use. If the average length of a rayrun is a lot smaller than the size of the partition, there will belittle gain, as the aggregate statistics associated with the partitionrepresent the entire partition and the number of ray runs that will berejected is reduced. Further, if the partition size is much larger thanthe ray run, the cost and complexity of checking for significant runs inthe partition is greater. Therefore, the cost of traversing thepartitioning must be balanced with the savings from finding voxels thatneed not be considered. If dealing with voxels that are relevant to theray is expensive, it is better to use a small partition size. If that isnot the case, a larger partition size can permit a faster traversalalgorithm.

A good compromise is to choose a run-length-partition size that is atleast as large as the short run length in the path of the ray, as wasdone when a partition size of 4 was chosen in FIG. 30. In this case, theray at each step will traverse at most two partitions. Therefore, theloop to traverse the run through the partitions can be unrolled into atwo-stage process, as is done in algorithm 3101.

Efficiency of Ray Traversal with a Partitioned Hierarchy of Runs

A three-dimensional run length partitioning of a volume provides anO(n³/p) compression of the volumetric data where p is the partition sizeand n³ is the size of the volume. To test a run against a partitioning,at most |(r _(j) −1)/p |+1 partition checks must be performed. The orderof the ray traversal algorithm is O(n/τ), where τ is the average runlength in the path of the ray. Therefore the ray traversal processthrough the partitioning is O(n/min(p, τ)). By choosing τ to be lessthan or equal to the partition size, the complexity of O(n/τ) ismaintained for the entire traversal process.

Conclusion

The foregoing Detailed Description has disclosed how to use theinventor's techniques for determining the voxels that are intersected bya ray using runs of cells that are intersected by projections of the rayon planes and for determining whether a ray intersects a voxel whoseassociated data will affect the ray and has further disclosed the bestmode known to the inventor of practicing his techniques. It will,however, be immediately apparent to those skilled in the relevanttechnologies that many other ways of practicing the techniques arepossible. For example, the techniques used to determine the cells of theplanes that are intersected by the projections are those described inthe parents of the present application; however, any other way ofdetermining the cells that are intersected by the projections could beused as well. Similarly, any technique of determining the intersectedvoxels from the intersected cells could be used in addition to thosedisclosed herein. With regard to the techniques for determining whethera ray run traverses a voxel of an encoding run that affects the ray, anymethod of encoding that distinguishes significant runs from theremainder of the encoding run may be used, and aggregate values may bemaintained for any useful partition of the encoding runs. Finally, anyuseful representation of the runs of cells, the runs of voxels, and theencoding runs may be employed.

Since that is the case, the Detailed Description is to be regarded asbeing in all respects exemplary and not restrictive, and the breadth ofthe invention disclosed herein is to be determined not from the DetailedDescription, but rather from the claims as interpreted with the fullbreadth permitted by the patent laws.

1. A method practiced in a computer system of determining voxels in anobject space that are intersected by a ray, the method comprising thesteps of: making projections of the ray on a plurality of planes in theobject space; determining cells in the planes that are intersected bythe projections; and using the intersected cells to determine theintersected voxels.
 2. The method set forth in claim 1 wherein the stepof determining cells in the planes includes the step of: for eachprojection determining a set of runs of cells that are intersected bythe projection.
 3. The method set forth in claim 2 wherein: the runs mayhave an order greater than
 1. 4. the method set forth in claim 2wherein: the object space has an axis that is the major axis relative tothe ray; and in the step of using the intersected cells, at thebeginning of the next run of voxels, if the first-order runs of cellsthat include the points in the projections corresponding to thebeginning of the next run of voxels end at the same major axiscoordinate, the first-order runs of cells together determine the nextrun of voxels through the ends of the first-order runs; and if the firstorder runs of cells that include the points in the projectionscorresponding to the beginning of the next run of voxels do not end atthe same major axis coordinate, the shorter first order run of cells andthe corresponding portion of the longer of the first order runs of cellsdetermine the next run of voxels, whereby the voxels intersected by theray are 26-connected.
 5. The method set forth in claim 4 wherein: in thestep of using the intersected cells, an extra cell is added to thebeginning of a first order run of cells prior to using the first orderrun of cells to determine a run of voxels, whereby the voxelsintersected by the ray are 6-connected.
 6. The method set forth in claim2 wherein: in the step of determining a set of runs, the set of runs fora given projection are determined in parallel.
 7. The method set forthin claim 2 wherein: the step of using the intersected cells to determinethe intersected voxels further includes the step of determining whetherthe intersected voxels are edge-connected or corner-connected.
 8. Themethod set forth in claim 7 wherein: in the step of determining whetherthe intersected voxels are edge-connected or corner connected, if one ofthe first-order runs of cells has a corner connection at a point and theother first order run of cells does not have a corner connection at thecorresponding point, the intersected voxels have an edge connection atthe corresponding point.
 9. The method set forth in claim 7 wherein: inthe step of determining whether the intersected voxels areedge-connected or corner connected, if both of the first-order runs ofcells have corner connections at a corresponding point, the intersectedvoxels have a corner connection at the corresponding point.
 10. Themethod set forth in claim 1 wherein: the object space has an axis thatis the major axis relative to the ray; and the plurality of planes istwo planes which intersect along the major axis.
 11. The method setforth in claim 10 wherein: the two planes intersect at right angles. 12.A method practiced in a computer system of traversing a volume with aparticular ray of a plurality thereof, the volume being subdivided intofirst runs of voxels, certain of the voxels being associated with datathat affects rays, and a ray intersecting one or more of the first runsand being defined as a set of second runs of voxels, and the methodcomprising the steps of: for a second run belonging to the particularray, determining whether the second run includes a voxel of a first runthat affects rays; and when the second run includes such a voxel,examining the associated data.
 13. The method set forth in claim 12wherein the first runs contain significant runs that include the certainvoxels; and the step of determining whether the particular ray's secondrun includes a voxel of a first run that affects rays includesdetermining whether the second run includes a voxel of a significantrun.
 14. The method set forth in claim 12 wherein: the volume has anaxis that is the major axis for both the particular ray and the firstruns of voxels.
 15. The method set forth in claim 14 wherein: there arethree sets of first runs, each set thereof having a different axis ofthe volume as its major axis.
 16. The method set forth in claim 12wherein: aggregate information is associated with partitions of thefirst runs, the aggregate information associated with a partitionindicating how one or more voxels in the partition affect rays; and inthe step of determining whether second run includes a voxel of a firstrun that affects rays, the aggregate information associated with apartition is used to determine whether the partition contains a voxelthat affects the particular ray.
 17. The method set forth in claim 16wherein: a first run has associated therewith a plurality of sets ofpartitions, the partitions in each set having a different length invoxels; and the step of determining includes the step of selecting oneof the sets of partitions in accordance with the lengths of the secondruns in the particular ray.
 18. The method set forth in claim 13wherein: aggregate information is associated with the significant runs,the aggregate information associated with the significant run indicatinghow one or more voxels in the partition affect rays; and the step ofdetermining whether second run includes a voxel of a first run thataffects rays includes using the aggregate information associated with asignificant run to determine whether the significant run contains avoxel that affects the particular ray.